These outputs will provide information on local areas, 2013 area units, and statistical area 2 units, and 2018 Meshblocks.
Note that all longitudes and latitudes should be in WGS1984 format
Four csv files:
#pragma nodebook off
#Use nodebook for better reproducibility https://github.com/uoa-eResearch/nodebook
%reload_ext nodebook.ipython
%nodebook disk phase2
# load libraries
import geopandas as gpd # vector data
import pandas as pd # tabular data, loading CSVs
import numpy as np # numeric data
from util import *
import matplotlib # plotting
from matplotlib_scalebar.scalebar import ScaleBar # scalebar for plot
import matplotlib.pyplot as plt # plotting
from tqdm.auto import tqdm # progress bars
tqdm.pandas()
import json
from shapely.geometry import Point, shape # creating points
import requests # web requests
plt.rcParams['figure.figsize'] = (20, 20)
ls()
| name | filesize (MB) | last modified | |
|---|---|---|---|
| 0 | 2013-mb-dataset-Total-New-Zealand-Household.csv | 37.12 | 2014-06-04 10:56:30.000000 |
| 1 | 2013-mb-dataset-Total-New-Zealand-Individual-Part-1.csv | 31.66 | 2014-04-01 10:13:15.197000 |
| 2 | 2018-census-electoral-population-meshblock-2020-data.csv | 5.60 | 2021-08-09 00:39:18.000000 |
| 3 | 2018_census_dwellings_by_SA2.xlsx | 0.43 | 2021-07-12 14:30:28.000000 |
| 4 | AC_Special_Housing_Area.zip | 0.29 | 2021-09-02 10:36:16.760000 |
| 5 | AucklandArea.gpkg | 0.09 | 2021-09-02 10:36:17.220000 |
| 6 | Geocoordinates_Direct_Transit_stops_AKL.xlsx | 0.01 | 2020-10-08 07:52:29.000000 |
| 7 | Individual_part1_totalNZ-wide_format_updated_16-7-20.csv | 36.58 | 2020-07-14 16:12:24.000000 |
| 8 | MASTER_UP_BaseZone_SHP.zip | 66.92 | 2021-07-19 02:23:51.137347 |
| 9 | Modified_Community_Boards_SHP.zip | 1.30 | 2021-07-19 02:16:07.650000 |
| 10 | area-unit-2013.gdb.zip | 14.21 | 2021-09-02 10:36:16.070000 |
| 11 | kx-nz-state-highway-on-ramps-off-ramps-SHP.zip | 0.14 | 2021-08-25 09:52:06.341291 |
| 12 | lds-nz-coastline-mean-high-water-FGDB.zip | 4.88 | 2021-08-31 16:25:16.250000 |
| 13 | lds-nz-coastlines-and-islands-polygons-topo-150k-FGDB.zip | 4.25 | 2021-08-05 15:57:49.020000 |
| 14 | lds-nz-primary-parcels-CLIPPED-4326.gpkg | 273.35 | 2021-09-02 10:37:32.410000 |
| 15 | lds-nz-primary-parcels-FGDB.zip | 79.03 | 2021-08-12 09:13:27.910000 |
| 16 | lds-nz-railway-centrelines-topo-150k-SHP.zip | 0.87 | 2021-08-31 10:01:46.315326 |
| 17 | lds-nz-road-centrelines-topo-150k-FGDB.zip | 36.06 | 2021-08-31 09:43:07.931704 |
| 18 | lds-nz-street-address-GPKG-CLIPPED.gpkg | 170.23 | 2021-09-02 10:37:19.490000 |
| 19 | lds-nz-street-address-GPKG.zip | 200.88 | 2021-09-02 10:37:26.070000 |
| 20 | meshblock-2013.gdb.zip | 106.57 | 2021-09-02 10:36:59.860000 |
| 21 | meshblock-2018-clipped-generalised.gdb.zip | 32.83 | 2021-09-02 10:36:26.010000 |
| 22 | nz-primary-parcels.gdb.zip | 55.00 | 2021-09-02 10:36:39.020000 |
| 23 | statsnzarea-unit-2013-FGDB.zip | 13.79 | 2021-08-31 16:56:17.150000 |
| 24 | statsnzmeshblock-higher-geographies-2018-generalised-FGDB.zip | 34.54 | 2021-08-05 14:53:54.350000 |
| 25 | statsnzpopulation-by-meshblock-2013-census-FGDB.zip | 82.11 | 2021-07-19 13:53:55.150631 |
| 26 | statsnzstatistical-area-2-higher-geographies-2018-clipped-generalis-FGDB.zip | 10.72 | 2021-08-30 17:12:55.591895 |
Total: 1299.0MB
df = gpd.read_file("input/Modified_Community_Boards_SHP.zip")
df.OBJECTID = df.OBJECTID.astype(int)
df = df.set_index("OBJECTID", drop=True)
df = df.sort_index()
display(df)
| Local_Area | geometry | |
|---|---|---|
| OBJECTID | ||
| 1 | Auckland Central | POLYGON ((1755802.315 5921956.091, 1755861.443... |
| 2 | Beach Haven-Birkenhead-Northcote | POLYGON ((1757287.966 5925962.738, 1757211.333... |
| 3 | Botany | POLYGON ((1770748.846 5912611.168, 1770889.893... |
| 4 | Devonport-Takapuna | POLYGON ((1755276.581 5932026.336, 1755278.305... |
| 5 | East Coast Bays | POLYGON ((1756125.006 5940268.048, 1756139.852... |
| 6 | Franklin-Beachlands-Hunua | MULTIPOLYGON (((1804302.354 5890738.079, 17905... |
| 7 | Franklin-Pukekohe | POLYGON ((1765085.620 5897344.807, 1765096.229... |
| 8 | Franklin-Waiuku | POLYGON ((1744829.308 5899882.633, 1744835.760... |
| 9 | Henderson-Massey | POLYGON ((1745963.138 5923457.510, 1745945.884... |
| 10 | Hibiscus Coast | MULTIPOLYGON (((1752023.352 5954803.281, 17520... |
| 11 | Howick | POLYGON ((1771447.460 5916900.636, 1771474.703... |
| 12 | Mangere-Otahuhu | MULTIPOLYGON (((1766142.746 5910838.845, 17661... |
| 13 | Manurewa | POLYGON ((1769432.488 5904664.673, 1769471.312... |
| 14 | Maungakiekie | POLYGON ((1759418.674 5915517.660, 1759596.351... |
| 15 | Mt Albert | POLYGON ((1756291.061 5918265.104, 1756266.530... |
| 16 | Mt Eden | POLYGON ((1751910.995 5920299.659, 1751917.632... |
| 17 | Orakei | POLYGON ((1759833.395 5920429.502, 1759837.143... |
| 18 | Otara | MULTIPOLYGON (((1765765.640 5909501.655, 17657... |
| 19 | Pakuranga | POLYGON ((1769523.742 5920466.582, 1769527.337... |
| 20 | Papakura | POLYGON ((1772188.699 5902202.011, 1772192.566... |
| 21 | Papatoetoe | POLYGON ((1765105.421 5908611.893, 1765113.661... |
| 22 | Puketapapa | POLYGON ((1753132.892 5915040.177, 1753141.101... |
| 23 | Rodney-Dairy Flat | POLYGON ((1746055.938 5948651.765, 1746067.112... |
| 24 | Rodney-Helensville | MULTIPOLYGON (((1710567.626 5967865.416, 17106... |
| 25 | Rodney-Kumeu-Riverhead | POLYGON ((1742532.808 5931237.574, 1742490.377... |
| 26 | Rodney-Warkworth | MULTIPOLYGON (((1761692.564 5984690.018, 17617... |
| 27 | Rodney-Wellsford | MULTIPOLYGON (((1746318.204 6000215.757, 17464... |
| 28 | Tamaki | POLYGON ((1765594.677 5917986.938, 1765609.808... |
| 29 | Titirangi | POLYGON ((1749038.038 5910572.842, 1749034.230... |
| 30 | Upper Harbour Local Board Area | MULTIPOLYGON (((1743880.554 5928979.708, 17438... |
| 31 | Waiheke | MULTIPOLYGON (((1793981.864 5931917.843, 17940... |
| 32 | Waitakere | MULTIPOLYGON (((1743890.310 5905057.188, 17438... |
| 33 | Whau | MULTIPOLYGON (((1748168.244 5916597.058, 17481... |
df = df.rename(columns={"Local_Area": "Name"})
df["Centroid_lon"] = df.centroid.to_crs(epsg=4326).x
df["Centroid_lat"] = df.centroid.to_crs(epsg=4326).y
df["Area"] = df.area
%%time
zones = gpd.read_file("input/MASTER_UP_BaseZone_SHP.zip")
zones
CPU times: user 12 s, sys: 626 ms, total: 12.6 s Wall time: 12.6 s
| OBJECTID | CONTOUR | created_da | DocumentUR | GlobalID | GROUPZONE | GROUPZONE_ | last_edite | NAME | PARCEL_BAS | ... | TYPE | TYPE_resol | VERSIONSTA | VERSIONS_1 | ZONE | ZONE_resol | ZONEHEIGHT | SHAPE_Leng | SHAPE_Area | geometry | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 1.0 | None | 20160718211 | None | {4C8F9436-7EA6-417E-B64F-15FCD44459F6} | 2 | Residential | 20161111010 | None | None | ... | None | None | 4 | Operative | 60 | Residential - Mixed Housing Urban Zone | NaN | 285.664016 | 2.050275e+03 | POLYGON ((1768030.306 5901206.846, 1768033.070... |
| 1 | 2.0 | None | 20160718211 | None | {604AAD87-8ED4-4111-8276-47CEE7E81F92} | 1 | Public Open Space | 20161111010 | None | None | ... | None | None | 4 | Operative | 33 | Open Space - Sport and Active Recreation Zone | NaN | 1246.837757 | 1.684599e+04 | POLYGON ((1764267.286 5919989.370, 1764218.153... |
| 2 | 3.0 | None | 20160718211 | None | {8D827DA8-BC5B-437A-B17A-532354F7D037} | 4 | Rural | 20161111010 | None | None | ... | None | None | 4 | Operative | 11 | Rural - Mixed Rural Zone | NaN | 3582.113246 | 6.841744e+05 | POLYGON ((1740091.195 5928308.839, 1740089.844... |
| 3 | 4.0 | None | 20160718211 | None | {96C9E266-3341-4C71-94F1-325F2EE45732} | 2 | Residential | 20161111010 | None | None | ... | None | None | 4 | Operative | 23 | Residential - Large Lot Zone | NaN | 317.098469 | 6.024226e+03 | POLYGON ((1750502.125 5928677.635, 1750456.131... |
| 4 | 5.0 | None | 20160718211 | None | {90B50FEE-45A3-4E88-819A-370751ACDE3D} | 1 | Public Open Space | 20161111010 | None | None | ... | None | None | 4 | Operative | 31 | Open Space - Conservation Zone | NaN | 230.836636 | 5.639794e+02 | POLYGON ((1741460.217 5918191.600, 1741476.745... |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 130295 | 130296.0 | None | 20161115151 | None | {B2F0FB45-80F6-41ED-AA57-A914B195B31E} | 1 | Public Open Space | 20161115151 | None | None | ... | None | None | 4 | Operative | 31 | Open Space - Conservation Zone | NaN | 179382.307787 | 9.040193e+07 | POLYGON ((1738072.631 5902593.421, 1738095.386... |
| 130296 | 130297.0 | None | 20161115151 | None | {A3F4EF61-9162-43C3-90BD-DC84D93C6A64} | 2 | Residential | 20161115151 | None | None | ... | None | None | 4 | Operative | 18 | Residential - Mixed Housing Suburban Zone | NaN | 2769.292040 | 2.920602e+05 | POLYGON ((1770189.642 5905789.775, 1770199.407... |
| 130297 | 130298.0 | None | 20161115151 | None | {A5FC7EB5-76E5-414C-96DD-715C671764B4} | 4 | Rural | 20161115151 | Kaipara South Head and Harbour coastal area | None | ... | None | None | 4 | Operative | 46 | Rural - Rural Coastal Zone | NaN | 27168.762393 | 6.344938e+06 | POLYGON ((1730253.401 5955767.484, 1730254.399... |
| 130298 | 130299.0 | None | 20161115151 | None | {C371B83C-35CE-46A1-B94F-F1E592F271F7} | 2 | Residential | 20161115151 | None | None | ... | None | None | 4 | Operative | 8 | Residential - Terrace Housing and Apartment Bu... | NaN | 1812.109533 | 1.536829e+05 | POLYGON ((1743032.759 5924130.564, 1742996.886... |
| 130299 | 130300.0 | None | 20161115151 | None | {31281924-C74D-4038-8260-FE6D886F4C94} | 2 | Residential | 20161115151 | None | None | ... | None | None | 4 | Operative | 18 | Residential - Mixed Housing Suburban Zone | NaN | 2511.267945 | 3.539606e+05 | POLYGON ((1748979.767 5925771.895, 1748986.615... |
130300 rows × 31 columns
zones[pd.isna(zones.ZONE_resol)]
| OBJECTID | CONTOUR | created_da | DocumentUR | GlobalID | GROUPZONE | GROUPZONE_ | last_edite | NAME | PARCEL_BAS | ... | TYPE | TYPE_resol | VERSIONSTA | VERSIONS_1 | ZONE | ZONE_resol | ZONEHEIGHT | SHAPE_Leng | SHAPE_Area | geometry | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 20853 | 20854.0 | None | 20160718211 | None | {2D83F680-9587-4D57-AF0B-CB7CA27C3D2C} | 6 | Special purpose zone | 20161111011 | None | None | ... | None | None | 4 | Operative | 58 | None | NaN | 72.629617 | 22.286945 | POLYGON ((1767684.860 5903714.023, 1767691.955... |
| 121765 | 121766.0 | None | 20160718211 | None | {FD3E8BB4-1979-42B5-9A77-5D04AE5190AF} | 6 | Special purpose zone | 20161111010 | None | None | ... | None | None | 4 | Operative | 58 | None | NaN | 86.135341 | 4.114962 | POLYGON ((1767684.083 5903713.690, 1767666.954... |
2 rows × 31 columns
zones.ZONE_resol = zones.ZONE_resol.fillna("Special")
zones.GROUPZONE_.value_counts()
General 52461 Coastal 42200 Residential 22398 Public Open Space 6667 Business 3237 Rural 2750 Special purpose zone 299 New growth 288 Name: GROUPZONE_, dtype: int64
zones.ZONE_resol.value_counts()
Road 47012 Coastal - General Coastal Marine Zone 26326 Coastal - Coastal Transition Zone 15703 Residential - Mixed Housing Suburban Zone 9775 Residential - Mixed Housing Urban Zone 5864 Strategic Transport Corridor Zone 4215 Residential - Single House Zone 3766 Open Space - Conservation Zone 3063 Open Space - Informal Recreation Zone 2926 Residential - Terrace Housing and Apartment Building Zone 2134 Rural - Rural Production Zone 1025 Water 1013 Business - Mixed Use Zone 782 Rural - Rural Coastal Zone 768 Business - Light Industry Zone 593 Business - Neighbourhood Centre Zone 564 Residential - Rural and Coastal Settlement Zone 480 Open Space - Sport and Active Recreation Zone 477 Business - Town Centre Zone 432 Residential - Large Lot Zone 379 Rural - Mixed Rural Zone 344 Rural - Countryside Living Zone 326 Future Urban Zone 282 Business - Local Centre Zone 266 Hauraki Gulf Islands 221 Business - City Centre Zone 191 Open Space - Community Zone 183 Rural - Waitakere Ranges Zone 172 Business - Metropolitan Centre Zone 160 Business - Heavy Industry Zone 130 Special Purpose - School Zone 126 Coastal - Mooring Zone 124 Business - General Business Zone 103 Rural - Rural Conservation Zone 71 Rural - Waitakere Foothills Zone 44 Special Purpose - Cemetery Zone 41 Special Purpose - Quarry Zone 34 Special Purpose - Maori Purpose Zone 32 Special Purpose - Major Recreation Facility Zone 28 Special Purpose - Healthcare Facility and Hospital Zone 23 Coastal - Marina Zone 20 Coastal - Ferry Terminal Zone 20 Open Space - Civic Spaces Zone 18 Business - Business Park Zone 16 Special Purpose - Airports and Airfields Zone 9 Green Infrastructure Corridor 6 Special Purpose - Tertiary Education Zone 4 Coastal - Minor Port Zone 4 Coastal - Defence Zone 3 Special 2 Name: ZONE_resol, dtype: int64
# Create a dictionary mapping the first n character of the zone description to the desired variable name
zones_of_interest = {
"Residential": "Residential_area",
"Residential - Single House Zone": "SH_area",
"Residential - Mixed Housing Suburban Zone": "MHS_area",
"Residential - Mixed Housing Urban Zone": "MHU_area",
"Residential - Terrace Housing and Apartment Building Zone": "THA_area",
"Residential - Large Lot Zone": "LL_area",
"Future Urban Zone": "FU_area",
"Hauraki Gulf Islands": "HGI_area",
"Business": "Business_area",
"Rural": "Rural_area",
"Open Space": "Open_area"
}
%%time
for k,v in tqdm(zones_of_interest.items()):
subset = zones[zones.ZONE_resol.str.startswith(k)] # Just the matched zones (residential, business, rural etc)
print(k, len(subset))
clipped = gpd.clip(df, subset) # Clip the local areas to the matched zones
df[v] = clipped.area # Store the resulting area. Unit is m²
Residential 22398 Residential - Single House Zone 3766 Residential - Mixed Housing Suburban Zone 9775 Residential - Mixed Housing Urban Zone 5864 Residential - Terrace Housing and Apartment Building Zone 2134 Residential - Large Lot Zone 379 Future Urban Zone 282 Hauraki Gulf Islands 221 Business 3237 Rural 2750 Open Space 6667 CPU times: user 5min 57s, sys: 247 ms, total: 5min 58s Wall time: 5min 57s
df = df.fillna(0)
skytower = Point(1757109.809, 5920500.841) # in NZGD2000 projection
df["Hdist_skytower"] = df.centroid.distance(skytower) # Unit is meters
df.plot(column="Hdist_skytower", legend=True)
<AxesSubplot:>
coastline = df.dissolve().boundary.iloc[0]
df["Coast_indicator"] = df.intersects(coastline)
df[["Name", "Coast_indicator"]]
| Name | Coast_indicator | |
|---|---|---|
| OBJECTID | ||
| 1 | Auckland Central | True |
| 2 | Beach Haven-Birkenhead-Northcote | True |
| 3 | Botany | True |
| 4 | Devonport-Takapuna | True |
| 5 | East Coast Bays | True |
| 6 | Franklin-Beachlands-Hunua | True |
| 7 | Franklin-Pukekohe | True |
| 8 | Franklin-Waiuku | True |
| 9 | Henderson-Massey | True |
| 10 | Hibiscus Coast | True |
| 11 | Howick | True |
| 12 | Mangere-Otahuhu | True |
| 13 | Manurewa | True |
| 14 | Maungakiekie | True |
| 15 | Mt Albert | False |
| 16 | Mt Eden | True |
| 17 | Orakei | True |
| 18 | Otara | True |
| 19 | Pakuranga | True |
| 20 | Papakura | True |
| 21 | Papatoetoe | True |
| 22 | Puketapapa | True |
| 23 | Rodney-Dairy Flat | True |
| 24 | Rodney-Helensville | True |
| 25 | Rodney-Kumeu-Riverhead | True |
| 26 | Rodney-Warkworth | True |
| 27 | Rodney-Wellsford | True |
| 28 | Tamaki | True |
| 29 | Titirangi | True |
| 30 | Upper Harbour Local Board Area | True |
| 31 | Waiheke | True |
| 32 | Waitakere | True |
| 33 | Whau | True |
households = pd.read_csv("input/2013-mb-dataset-Total-New-Zealand-Household.csv", encoding='unicode_escape')
households
/home/nyou045/git/nodebook/nodebook/nodebookcore.py:283: DtypeWarning: Columns (1,2) have mixed types.Specify dtype option on import or set low_memory=False. exec(compile(block, '<string>', mode='exec'), env)
| Area_Code_and_Description | Code | Description | 2001_Census_total_households_in_occupied_private_dwellings | 2006_Census_total_households_in_occupied_private_dwellings | 2013_Census_total_households_in_occupied_private_dwellings | 2001_Census_household_composition_for_households_in_occupied_private_dwellings_One-family_household_(with_or_without_other_people) | 2001_Census_household_composition_for_households_in_occupied_private_dwellings_Two-family_household_(with_or_without_other_people) | 2001_Census_household_composition_for_households_in_occupied_private_dwellings_Three_or_more_family_household_(with_or_without_other_people) | 2001_Census_household_composition_for_households_in_occupied_private_dwellings_Other_multi-person_household | ... | 2006_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Not_Elsewhere_Included(22) | 2006_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Total_households_in_occupied_private_dwellings | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_No_Access_to_Telecommunication_Systems | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Access_to_a_Cellphone/Mobile_Phone | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Access_to_a_Telephone | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Access_to_a_Fax_Machine | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Access_to_the_Internet | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Total_households_stated | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Not_Elsewhere_Included(22) | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Total_households_in_occupied_private_dwellings | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | MB 0000100 | 100 | NaN | 3.0 | 3.0 | 3.0 | ..C | ..C | ..C | ..C | ... | ..C | 3.0 | ..C | ..C | ..C | ..C | ..C | ..C | ..C | 3.0 |
| 1 | MB 0000200 | 200 | NaN | 27.0 | 27.0 | 30.0 | 15 | 3 | 0 | 0 | ... | 3 | 24.0 | 3 | 9 | 24 | 6 | 15 | 30 | 3 | 30.0 |
| 2 | MB 0000300 | 300 | NaN | 21.0 | 24.0 | 21.0 | 18 | 3 | 0 | 0 | ... | 0 | 24.0 | 0 | 6 | 12 | 3 | 9 | 18 | 3 | 18.0 |
| 3 | MB 0000400 | 400 | NaN | 9.0 | 9.0 | 9.0 | 6 | ..C | ..C | ..C | ... | ..C | 12.0 | ..C | ..C | 9 | ..C | 6 | 9 | ..C | 9.0 |
| 4 | MB 0000501 | 501 | NaN | 0.0 | 0.0 | 0.0 | ..C | ..C | ..C | ..C | ... | ..C | 0.0 | ..C | ..C | ..C | ..C | ..C | ..C | ..C | 0.0 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 49038 | Footnotes for specific variables or categories | 23 | Median total household income is rounded to th... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
| 49039 | Symbols | ..C | Confidential | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
| 49040 | Symbols | .. | not available | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
| 49041 | Symbols | * | not able to be calculated | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
| 49042 | Source | Statistics New Zealand | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
49043 rows × 258 columns
key = "2013_Census_total_household_income_(grouped)(2)(3)(4)_for_households_in_occupied_private_dwellings_Median_household_income_($)(18)(23)"
households_with_income = households[households.Area_Code_and_Description.str.startswith("MB") & (households[key] != "..C") & (households[key] != "*")].copy()
print(len(households_with_income))
households_with_income[key] = households_with_income[key].astype(int)
households_with_income[key].plot(kind="hist", bins=20, figsize=(10,10))
40300
<AxesSubplot:ylabel='Frequency'>
meshblocks = gpd.read_file("input/statsnzpopulation-by-meshblock-2013-census-FGDB.zip!population-by-meshblock-2013-census.gdb")
meshblocks
| Meshblock | MeshblockNumber | Population_Count_Usual_Resident_2013 | Population_Count_Census_Night_2013 | geometry | |
|---|---|---|---|---|---|
| 0 | MB 0352700 | 0352700 | 0 | 0 | MULTIPOLYGON (((1753237.550 5923918.395, 17532... |
| 1 | MB 0728500 | 0728500 | 9 | 6 | MULTIPOLYGON (((1761011.438 5905840.848, 17610... |
| 2 | MB 0829300 | 0829300 | 93 | 99 | MULTIPOLYGON (((1739739.656 5899165.556, 17397... |
| 3 | MB 1280801 | 1280801 | 0 | 0 | MULTIPOLYGON (((1869852.595 5695096.558, 18714... |
| 4 | MB 2360001 | 2360001 | 0 | 0 | MULTIPOLYGON (((1623601.773 5423210.098, 16235... |
| ... | ... | ... | ... | ... | ... |
| 46616 | MB 0074002 | 0074002 | 0 | 0 | MULTIPOLYGON (((1741178.600 6046572.236, 17414... |
| 46617 | MB 3208002 | 3208002 | 78 | 81 | MULTIPOLYGON (((1770892.430 5911519.906, 17708... |
| 46618 | MB 3208003 | 3208003 | 36 | 36 | MULTIPOLYGON (((1771025.156 5911674.629, 17709... |
| 46619 | MB 3209001 | 3209001 | 84 | 87 | MULTIPOLYGON (((1771731.425 5912665.799, 17717... |
| 46620 | MB 3209002 | 3209002 | 75 | 81 | MULTIPOLYGON (((1771758.103 5912422.760, 17717... |
46621 rows × 5 columns
meshblocks.geometry = meshblocks.representative_point()
meshblocks = meshblocks.merge(households, left_on="Meshblock", right_on="Area_Code_and_Description")
meshblocks = gpd.sjoin(df, meshblocks, op="intersects")
meshblocks
| Name | geometry | Centroid_lon | Centroid_lat | Area | Residential_area | SH_area | MHS_area | MHU_area | THA_area | ... | 2006_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Not_Elsewhere_Included(22) | 2006_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Total_households_in_occupied_private_dwellings | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_No_Access_to_Telecommunication_Systems | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Access_to_a_Cellphone/Mobile_Phone | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Access_to_a_Telephone | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Access_to_a_Fax_Machine | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Access_to_the_Internet | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Total_households_stated | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Not_Elsewhere_Included(22) | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Total_households_in_occupied_private_dwellings | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| OBJECTID | |||||||||||||||||||||
| 1 | Auckland Central | POLYGON ((1755802.315 5921956.091, 1755861.443... | 174.753763 | -36.855194 | 1.942334e+07 | 6.114247e+06 | 3.264382e+06 | 1.005111e+06 | 9.761411e+05 | 8.686144e+05 | ... | ..C | 0.0 | ..C | ..C | ..C | ..C | ..C | ..C | ..C | 0.0 |
| 1 | Auckland Central | POLYGON ((1755802.315 5921956.091, 1755861.443... | 174.753763 | -36.855194 | 1.942334e+07 | 6.114247e+06 | 3.264382e+06 | 1.005111e+06 | 9.761411e+05 | 8.686144e+05 | ... | ..C | 0.0 | ..C | ..C | ..C | ..C | ..C | ..C | ..C | 0.0 |
| 1 | Auckland Central | POLYGON ((1755802.315 5921956.091, 1755861.443... | 174.753763 | -36.855194 | 1.942334e+07 | 6.114247e+06 | 3.264382e+06 | 1.005111e+06 | 9.761411e+05 | 8.686144e+05 | ... | ..C | 0.0 | ..C | ..C | ..C | ..C | ..C | ..C | ..C | 0.0 |
| 1 | Auckland Central | POLYGON ((1755802.315 5921956.091, 1755861.443... | 174.753763 | -36.855194 | 1.942334e+07 | 6.114247e+06 | 3.264382e+06 | 1.005111e+06 | 9.761411e+05 | 8.686144e+05 | ... | ..C | 6.0 | 0 | 15 | 6 | 0 | 18 | 18 | 6 | 24.0 |
| 1 | Auckland Central | POLYGON ((1755802.315 5921956.091, 1755861.443... | 174.753763 | -36.855194 | 1.942334e+07 | 6.114247e+06 | 3.264382e+06 | 1.005111e+06 | 9.761411e+05 | 8.686144e+05 | ... | ..C | 9.0 | ..C | 6 | 6 | ..C | 9 | 6 | ..C | 9.0 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 33 | Whau | MULTIPOLYGON (((1748168.244 5916597.058, 17481... | 174.684472 | -36.906902 | 2.682152e+07 | 1.602986e+07 | 9.520285e+05 | 6.562767e+06 | 6.248986e+06 | 2.009848e+06 | ... | ..C | 0.0 | ..C | ..C | ..C | ..C | ..C | ..C | ..C | 0.0 |
| 33 | Whau | MULTIPOLYGON (((1748168.244 5916597.058, 17481... | 174.684472 | -36.906902 | 2.682152e+07 | 1.602986e+07 | 9.520285e+05 | 6.562767e+06 | 6.248986e+06 | 2.009848e+06 | ... | 3 | 42.0 | 3 | 39 | 33 | 0 | 30 | 45 | 0 | 45.0 |
| 33 | Whau | MULTIPOLYGON (((1748168.244 5916597.058, 17481... | 174.684472 | -36.906902 | 2.682152e+07 | 1.602986e+07 | 9.520285e+05 | 6.562767e+06 | 6.248986e+06 | 2.009848e+06 | ... | 3 | 27.0 | 0 | 9 | 12 | 3 | 9 | 12 | 6 | 18.0 |
| 33 | Whau | MULTIPOLYGON (((1748168.244 5916597.058, 17481... | 174.684472 | -36.906902 | 2.682152e+07 | 1.602986e+07 | 9.520285e+05 | 6.562767e+06 | 6.248986e+06 | 2.009848e+06 | ... | ..C | 3.0 | ..C | ..C | ..C | ..C | ..C | ..C | ..C | 3.0 |
| 33 | Whau | MULTIPOLYGON (((1748168.244 5916597.058, 17481... | 174.684472 | -36.906902 | 2.682152e+07 | 1.602986e+07 | 9.520285e+05 | 6.562767e+06 | 6.248986e+06 | 2.009848e+06 | ... | ..C | 0.0 | ..C | ..C | ..C | ..C | ..C | ..C | ..C | 0.0 |
11527 rows × 281 columns
group = meshblocks.groupby("OBJECTID")
df["Census2013_population"] = group["Population_Count_Usual_Resident_2013"].sum()
df["Census2013_dwellings"] = group["2013_Census_total_households_in_occupied_private_dwellings"].sum()
meshblocks_with_income = meshblocks[(meshblocks[key] != "..C") & (meshblocks[key] != "*")].copy()
meshblocks_with_income[key] = meshblocks_with_income[key].astype(int)
meshblocks_with_income
| Name | geometry | Centroid_lon | Centroid_lat | Area | Residential_area | SH_area | MHS_area | MHU_area | THA_area | ... | 2006_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Not_Elsewhere_Included(22) | 2006_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Total_households_in_occupied_private_dwellings | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_No_Access_to_Telecommunication_Systems | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Access_to_a_Cellphone/Mobile_Phone | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Access_to_a_Telephone | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Access_to_a_Fax_Machine | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Access_to_the_Internet | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Total_households_stated | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Not_Elsewhere_Included(22) | 2013_Census_access_to_telecommunications(20)(21)_for_households_in_occupied_private_dwellings_Total_households_in_occupied_private_dwellings | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| OBJECTID | |||||||||||||||||||||
| 1 | Auckland Central | POLYGON ((1755802.315 5921956.091, 1755861.443... | 174.753763 | -36.855194 | 1.942334e+07 | 6.114247e+06 | 3.264382e+06 | 1.005111e+06 | 9.761411e+05 | 8.686144e+05 | ... | ..C | 6.0 | 0 | 15 | 6 | 0 | 18 | 18 | 6 | 24.0 |
| 1 | Auckland Central | POLYGON ((1755802.315 5921956.091, 1755861.443... | 174.753763 | -36.855194 | 1.942334e+07 | 6.114247e+06 | 3.264382e+06 | 1.005111e+06 | 9.761411e+05 | 8.686144e+05 | ... | ..C | 9.0 | ..C | 6 | 6 | ..C | 9 | 6 | ..C | 9.0 |
| 1 | Auckland Central | POLYGON ((1755802.315 5921956.091, 1755861.443... | 174.753763 | -36.855194 | 1.942334e+07 | 6.114247e+06 | 3.264382e+06 | 1.005111e+06 | 9.761411e+05 | 8.686144e+05 | ... | 3 | 51.0 | 0 | 45 | 36 | 6 | 45 | 51 | 0 | 51.0 |
| 1 | Auckland Central | POLYGON ((1755802.315 5921956.091, 1755861.443... | 174.753763 | -36.855194 | 1.942334e+07 | 6.114247e+06 | 3.264382e+06 | 1.005111e+06 | 9.761411e+05 | 8.686144e+05 | ... | 0 | 48.0 | 3 | 45 | 36 | 6 | 39 | 48 | 0 | 48.0 |
| 1 | Auckland Central | POLYGON ((1755802.315 5921956.091, 1755861.443... | 174.753763 | -36.855194 | 1.942334e+07 | 6.114247e+06 | 3.264382e+06 | 1.005111e+06 | 9.761411e+05 | 8.686144e+05 | ... | 0 | 36.0 | 0 | 27 | 24 | 6 | 27 | 27 | 0 | 30.0 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 33 | Whau | MULTIPOLYGON (((1748168.244 5916597.058, 17481... | 174.684472 | -36.906902 | 2.682152e+07 | 1.602986e+07 | 9.520285e+05 | 6.562767e+06 | 6.248986e+06 | 2.009848e+06 | ... | 0 | 54.0 | 3 | 39 | 45 | 6 | 42 | 51 | 3 | 54.0 |
| 33 | Whau | MULTIPOLYGON (((1748168.244 5916597.058, 17481... | 174.684472 | -36.906902 | 2.682152e+07 | 1.602986e+07 | 9.520285e+05 | 6.562767e+06 | 6.248986e+06 | 2.009848e+06 | ... | 6 | 66.0 | 3 | 48 | 51 | 6 | 45 | 66 | 6 | 69.0 |
| 33 | Whau | MULTIPOLYGON (((1748168.244 5916597.058, 17481... | 174.684472 | -36.906902 | 2.682152e+07 | 1.602986e+07 | 9.520285e+05 | 6.562767e+06 | 6.248986e+06 | 2.009848e+06 | ... | 3 | 54.0 | 0 | 48 | 42 | 3 | 45 | 54 | 6 | 60.0 |
| 33 | Whau | MULTIPOLYGON (((1748168.244 5916597.058, 17481... | 174.684472 | -36.906902 | 2.682152e+07 | 1.602986e+07 | 9.520285e+05 | 6.562767e+06 | 6.248986e+06 | 2.009848e+06 | ... | 3 | 42.0 | 3 | 39 | 33 | 0 | 30 | 45 | 0 | 45.0 |
| 33 | Whau | MULTIPOLYGON (((1748168.244 5916597.058, 17481... | 174.684472 | -36.906902 | 2.682152e+07 | 1.602986e+07 | 9.520285e+05 | 6.562767e+06 | 6.248986e+06 | 2.009848e+06 | ... | 3 | 27.0 | 0 | 9 | 12 | 3 | 9 | 12 | 6 | 18.0 |
10830 rows × 281 columns
meshblocks_with_income["income*dwellings"] = meshblocks_with_income[key] * meshblocks_with_income["2013_Census_total_households_in_occupied_private_dwellings"]
meshblocks_with_income["income*dwellings"]
OBJECTID
1 1240800.0
1 1350000.0
1 4717500.0
1 4440000.0
1 2550000.0
...
33 2565000.0
33 3622500.0
33 3900000.0
33 1836000.0
33 930600.0
Name: income*dwellings, Length: 10830, dtype: float64
group = meshblocks_with_income.groupby("OBJECTID")
df["Census2013_avg_HH_income"] = group["income*dwellings"].sum() / group["2013_Census_total_households_in_occupied_private_dwellings"].sum()
df[["Name", "Census2013_population", "Census2013_dwellings", "Census2013_avg_HH_income"]].sort_values(by="Census2013_population", ascending=False)
| Name | Census2013_population | Census2013_dwellings | Census2013_avg_HH_income | |
|---|---|---|---|---|
| OBJECTID | ||||
| 9 | Henderson-Massey | 107658 | 34461.0 | 68347.842761 |
| 2 | Beach Haven-Birkenhead-Northcote | 82431 | 28344.0 | 80154.552207 |
| 13 | Manurewa | 82230 | 22584.0 | 69202.861326 |
| 17 | Orakei | 79581 | 29001.0 | 110143.425140 |
| 1 | Auckland Central | 77058 | 31506.0 | 86499.398108 |
| 33 | Whau | 72585 | 23928.0 | 65948.529781 |
| 12 | Mangere-Otahuhu | 70998 | 17427.0 | 61855.553636 |
| 4 | Devonport-Takapuna | 55398 | 20307.0 | 89560.517751 |
| 30 | Upper Harbour Local Board Area | 53571 | 17079.0 | 92354.882984 |
| 22 | Puketapapa | 52959 | 16695.0 | 75460.348858 |
| 15 | Mt Albert | 47913 | 15630.0 | 97670.686671 |
| 3 | Botany | 46878 | 13893.0 | 93325.464764 |
| 16 | Mt Eden | 46752 | 16155.0 | 83313.629602 |
| 20 | Papakura | 45687 | 14919.0 | 68043.070707 |
| 21 | Papatoetoe | 45633 | 13272.0 | 65055.631090 |
| 10 | Hibiscus Coast | 45141 | 17397.0 | 70786.765722 |
| 5 | East Coast Bays | 44634 | 15354.0 | 96602.171362 |
| 28 | Tamaki | 42189 | 13551.0 | 64390.277469 |
| 11 | Howick | 40299 | 13836.0 | 85281.267636 |
| 19 | Pakuranga | 39924 | 13206.0 | 82096.046353 |
| 29 | Titirangi | 35958 | 12192.0 | 80713.729357 |
| 7 | Franklin-Pukekohe | 31176 | 10863.0 | 78094.119279 |
| 18 | Otara | 30033 | 6693.0 | 58229.733394 |
| 14 | Maungakiekie | 27780 | 10368.0 | 76915.212333 |
| 6 | Franklin-Beachlands-Hunua | 21111 | 7332.0 | 97416.701031 |
| 24 | Rodney-Helensville | 17832 | 6381.0 | 78202.750119 |
| 26 | Rodney-Warkworth | 17637 | 6912.0 | 62686.367596 |
| 8 | Franklin-Waiuku | 13581 | 5046.0 | 68150.358852 |
| 32 | Waitakere | 12477 | 4449.0 | 84773.121192 |
| 31 | Waiheke | 8310 | 3597.0 | 55066.946309 |
| 25 | Rodney-Kumeu-Riverhead | 7257 | 2472.0 | 86943.048780 |
| 23 | Rodney-Dairy Flat | 6564 | 2130.0 | 106361.971831 |
| 27 | Rodney-Wellsford | 5595 | 2115.0 | 51240.028490 |
df
| Name | geometry | Centroid_lon | Centroid_lat | Area | Residential_area | SH_area | MHS_area | MHU_area | THA_area | ... | FU_area | HGI_area | Business_area | Rural_area | Open_area | Hdist_skytower | Coast_indicator | Census2013_population | Census2013_dwellings | Census2013_avg_HH_income | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| OBJECTID | |||||||||||||||||||||
| 1 | Auckland Central | POLYGON ((1755802.315 5921956.091, 1755861.443... | 174.753763 | -36.855194 | 1.942334e+07 | 6.114247e+06 | 3.264382e+06 | 1.005111e+06 | 9.761411e+05 | 8.686144e+05 | ... | 0.000000e+00 | 0.000000e+00 | 4.770257e+06 | 0.000000e+00 | 2.505081e+06 | 1067.517983 | True | 77058 | 31506.0 | 86499.398108 |
| 2 | Beach Haven-Birkenhead-Northcote | POLYGON ((1757287.966 5925962.738, 1757211.333... | 174.721897 | -36.794764 | 3.406603e+07 | 1.989760e+07 | 5.073543e+06 | 1.007797e+07 | 3.677952e+06 | 9.609651e+05 | ... | 0.000000e+00 | 0.000000e+00 | 2.554498e+06 | 0.000000e+00 | 6.442079e+06 | 6948.094501 | True | 82431 | 28344.0 | 80154.552207 |
| 3 | Botany | POLYGON ((1770748.846 5912611.168, 1770889.893... | 174.917904 | -36.957505 | 3.978830e+07 | 1.437879e+07 | 4.003298e+05 | 9.023214e+06 | 2.665899e+06 | 1.402369e+06 | ... | 0.000000e+00 | 0.000000e+00 | 6.222613e+06 | 9.837755e+06 | 4.056828e+06 | 18420.075506 | True | 46878 | 13893.0 | 93325.464764 |
| 4 | Devonport-Takapuna | POLYGON ((1755276.581 5932026.336, 1755278.305... | 174.769923 | -36.789009 | 2.111121e+07 | 1.240583e+07 | 2.595033e+06 | 6.760630e+06 | 2.528557e+06 | 5.216074e+05 | ... | 0.000000e+00 | 0.000000e+00 | 1.165067e+06 | 0.000000e+00 | 2.289014e+06 | 6621.239572 | True | 55398 | 20307.0 | 89560.517751 |
| 5 | East Coast Bays | POLYGON ((1756125.006 5940268.048, 1756139.852... | 174.733122 | -36.701804 | 3.038568e+07 | 1.402263e+07 | 2.429714e+06 | 8.490284e+06 | 9.649423e+05 | 2.219014e+05 | ... | 0.000000e+00 | 0.000000e+00 | 3.126271e+05 | 7.860819e+06 | 3.829094e+06 | 16467.095662 | True | 44634 | 15354.0 | 96602.171362 |
| 6 | Franklin-Beachlands-Hunua | MULTIPOLYGON (((1804302.354 5890738.079, 17905... | 175.101754 | -37.043410 | 7.778987e+08 | 7.835879e+06 | 6.069391e+06 | 3.296486e+04 | 0.000000e+00 | 3.034730e+04 | ... | 1.048487e+07 | 0.000000e+00 | 3.582480e+06 | 4.948783e+08 | 8.000078e+07 | 37193.822526 | True | 21111 | 7332.0 | 97416.701031 |
| 7 | Franklin-Pukekohe | POLYGON ((1765085.620 5897344.807, 1765096.229... | 174.864884 | -37.151426 | 2.717220e+08 | 1.622679e+07 | 6.467033e+06 | 5.660850e+06 | 3.079772e+06 | 3.792095e+05 | ... | 2.256611e+07 | 0.000000e+00 | 1.981925e+06 | 2.111716e+08 | 4.301604e+06 | 34851.656145 | True | 31176 | 10863.0 | 78094.119279 |
| 8 | Franklin-Waiuku | POLYGON ((1744829.308 5899882.633, 1744835.760... | 174.666624 | -37.164727 | 2.997611e+08 | 7.599913e+06 | 4.391057e+05 | 2.379868e+06 | 1.111391e+05 | 0.000000e+00 | ... | 1.885628e+05 | 0.000000e+00 | 5.024976e+06 | 2.704684e+08 | 4.015675e+06 | 36121.662648 | True | 13581 | 5046.0 | 68150.358852 |
| 9 | Henderson-Massey | POLYGON ((1745963.138 5923457.510, 1745945.884... | 174.621864 | -36.853375 | 5.321961e+07 | 3.041467e+07 | 5.403257e+06 | 6.196787e+06 | 1.438878e+07 | 4.167398e+06 | ... | 2.356178e+06 | 0.000000e+00 | 5.410361e+06 | 6.572595e+05 | 5.824062e+06 | 12525.129414 | True | 107658 | 34461.0 | 68347.842761 |
| 10 | Hibiscus Coast | MULTIPOLYGON (((1752023.352 5954803.281, 17520... | 174.720160 | -36.618022 | 7.969308e+07 | 2.547122e+07 | 1.531569e+07 | 1.601068e+06 | 1.640568e+06 | 4.426904e+05 | ... | 3.334191e+06 | 0.000000e+00 | 3.122063e+06 | 2.421639e+07 | 1.334561e+07 | 25832.060623 | True | 45141 | 17397.0 | 70786.765722 |
| 11 | Howick | POLYGON ((1771447.460 5916900.636, 1771474.703... | 174.926939 | -36.902134 | 1.461575e+07 | 1.060461e+07 | 3.809348e+06 | 5.160719e+06 | 1.340626e+06 | 2.939136e+05 | ... | 0.000000e+00 | 0.000000e+00 | 3.478985e+05 | 5.208255e+04 | 1.331571e+06 | 15854.353763 | True | 40299 | 13836.0 | 85281.267636 |
| 12 | Mangere-Otahuhu | MULTIPOLYGON (((1766142.746 5910838.845, 17661... | 174.794437 | -36.972479 | 5.250979e+07 | 1.304669e+07 | 1.231326e+06 | 6.294665e+06 | 3.943621e+06 | 1.366092e+06 | ... | 1.146902e+06 | 0.000000e+00 | 8.392959e+06 | 2.711786e+06 | 7.736789e+06 | 14070.305771 | True | 70998 | 17427.0 | 61855.553636 |
| 13 | Manurewa | POLYGON ((1769432.488 5904664.673, 1769471.312... | 174.886011 | -37.018487 | 3.711616e+07 | 1.691256e+07 | 1.083351e+06 | 1.234376e+07 | 2.890852e+06 | 5.945951e+05 | ... | 1.208450e+06 | 0.000000e+00 | 6.011804e+06 | 3.124826e+05 | 6.021889e+06 | 21865.784349 | True | 82230 | 22584.0 | 69202.861326 |
| 14 | Maungakiekie | POLYGON ((1759418.674 5915517.660, 1759596.351... | 174.799976 | -36.915243 | 1.740899e+07 | 5.002987e+06 | 3.103315e+05 | 2.479102e+06 | 1.117549e+06 | 1.096004e+06 | ... | 0.000000e+00 | 0.000000e+00 | 6.424025e+06 | 0.000000e+00 | 2.818621e+06 | 8151.195589 | True | 27780 | 10368.0 | 76915.212333 |
| 15 | Mt Albert | POLYGON ((1756291.061 5918265.104, 1756266.530... | 174.765299 | -36.884922 | 1.385353e+07 | 8.308876e+06 | 3.715019e+06 | 2.058002e+06 | 1.566702e+06 | 9.691525e+05 | ... | 0.000000e+00 | 0.000000e+00 | 1.219666e+06 | 0.000000e+00 | 9.752424e+05 | 4067.283622 | False | 47913 | 15630.0 | 97670.686671 |
| 16 | Mt Eden | POLYGON ((1751910.995 5920299.659, 1751917.632... | 174.718360 | -36.879882 | 1.448821e+07 | 8.270683e+06 | 1.092566e+06 | 2.464966e+06 | 3.332548e+06 | 1.380604e+06 | ... | 0.000000e+00 | 0.000000e+00 | 1.026338e+06 | 0.000000e+00 | 1.709990e+06 | 5244.560531 | True | 46752 | 16155.0 | 83313.629602 |
| 17 | Orakei | POLYGON ((1759833.395 5920429.502, 1759837.143... | 174.830542 | -36.870391 | 3.234712e+07 | 1.873652e+07 | 1.811262e+06 | 1.210010e+07 | 3.663284e+06 | 1.161878e+06 | ... | 0.000000e+00 | 0.000000e+00 | 1.326911e+06 | 0.000000e+00 | 5.326449e+06 | 6567.770945 | True | 79581 | 29001.0 | 110143.425140 |
| 18 | Otara | MULTIPOLYGON (((1765765.640 5909501.655, 17657... | 174.883562 | -36.963798 | 1.213648e+07 | 5.101724e+06 | 1.985154e+04 | 1.887588e+06 | 2.723733e+06 | 4.705511e+05 | ... | 0.000000e+00 | 0.000000e+00 | 2.434656e+06 | 0.000000e+00 | 2.367666e+06 | 16766.406523 | True | 30033 | 6693.0 | 58229.733394 |
| 19 | Pakuranga | POLYGON ((1769523.742 5920466.582, 1769527.337... | 174.891901 | -36.898057 | 1.527965e+07 | 9.096320e+06 | 3.916487e+04 | 7.010512e+06 | 1.552850e+06 | 4.937933e+05 | ... | 0.000000e+00 | 0.000000e+00 | 5.185541e+05 | 0.000000e+00 | 2.768248e+06 | 12812.703793 | True | 39924 | 13206.0 | 82096.046353 |
| 20 | Papakura | POLYGON ((1772188.699 5902202.011, 1772192.566... | 174.940390 | -37.059940 | 4.022528e+07 | 1.808030e+07 | 1.423108e+06 | 1.187048e+07 | 4.001187e+06 | 2.503432e+05 | ... | 6.270214e+06 | 0.000000e+00 | 4.066388e+06 | 2.638487e+06 | 2.646971e+06 | 28340.533206 | True | 45687 | 14919.0 | 68043.070707 |
| 21 | Papatoetoe | POLYGON ((1765105.421 5908611.893, 1765113.661... | 174.849291 | -36.987020 | 2.497254e+07 | 9.380759e+06 | 1.026630e+06 | 6.202250e+06 | 1.652101e+06 | 4.997783e+05 | ... | 6.587091e+05 | 0.000000e+00 | 7.079133e+06 | 0.000000e+00 | 3.473963e+06 | 17234.007915 | True | 45633 | 13272.0 | 65055.631090 |
| 22 | Puketapapa | POLYGON ((1753132.892 5915040.177, 1753141.101... | 174.741502 | -36.916029 | 1.871805e+07 | 1.108573e+07 | 7.062543e+05 | 6.066849e+06 | 2.942080e+06 | 1.360481e+06 | ... | 0.000000e+00 | 0.000000e+00 | 1.061183e+06 | 0.000000e+00 | 3.279811e+06 | 7732.493593 | True | 52959 | 16695.0 | 75460.348858 |
| 23 | Rodney-Dairy Flat | POLYGON ((1746055.938 5948651.765, 1746067.112... | 174.636168 | -36.672259 | 1.132876e+08 | 2.928810e+06 | 4.415845e+05 | 1.082742e+06 | 7.851764e+05 | 1.177736e+05 | ... | 2.571329e+07 | 0.000000e+00 | 2.402994e+05 | 7.596025e+07 | 3.594364e+06 | 22546.984200 | True | 6564 | 2130.0 | 106361.971831 |
| 24 | Rodney-Helensville | MULTIPOLYGON (((1710567.626 5967865.416, 17106... | 174.425342 | -36.646179 | 7.974650e+08 | 5.122441e+06 | 2.915926e+06 | 4.710662e+05 | 8.037448e+04 | 0.000000e+00 | ... | 4.384484e+05 | 0.000000e+00 | 5.046133e+05 | 7.218524e+08 | 2.634552e+07 | 37521.293694 | True | 17832 | 6381.0 | 78202.750119 |
| 25 | Rodney-Kumeu-Riverhead | POLYGON ((1742532.808 5931237.574, 1742490.377... | 174.538408 | -36.780645 | 5.536146e+07 | 4.569329e+06 | 3.637463e+06 | 8.916239e+05 | 4.024164e+04 | 0.000000e+00 | ... | 8.770729e+06 | 0.000000e+00 | 1.003907e+06 | 3.686394e+07 | 8.429452e+05 | 21330.923764 | True | 7257 | 2472.0 | 86943.048780 |
| 26 | Rodney-Warkworth | MULTIPOLYGON (((1761692.564 5984690.018, 17617... | 174.613063 | -36.444991 | 6.684679e+08 | 1.681370e+07 | 5.952212e+06 | 1.799798e+05 | 6.208502e+04 | 0.000000e+00 | ... | 1.266006e+07 | 0.000000e+00 | 1.633665e+06 | 5.541282e+08 | 3.901882e+07 | 46699.922300 | True | 17637 | 6912.0 | 62686.367596 |
| 27 | Rodney-Wellsford | MULTIPOLYGON (((1746318.204 6000215.757, 17464... | 174.538530 | -36.291690 | 6.404974e+08 | 1.680348e+06 | 1.152092e+06 | 0.000000e+00 | 0.000000e+00 | 0.000000e+00 | ... | 1.064854e+06 | 0.000000e+00 | 4.148911e+05 | 5.748374e+08 | 3.130978e+07 | 64927.125388 | True | 5595 | 2115.0 | 51240.028490 |
| 28 | Tamaki | POLYGON ((1765594.677 5917986.938, 1765609.808... | 174.848247 | -36.901670 | 1.896657e+07 | 8.372005e+06 | 1.763932e+04 | 3.137190e+06 | 2.703557e+06 | 2.513619e+06 | ... | 0.000000e+00 | 0.000000e+00 | 5.024667e+06 | 0.000000e+00 | 2.233061e+06 | 9688.801624 | True | 42189 | 13551.0 | 64390.277469 |
| 29 | Titirangi | POLYGON ((1749038.038 5910572.842, 1749034.230... | 174.629790 | -36.939520 | 4.058868e+07 | 1.476243e+07 | 1.672204e+06 | 2.303941e+06 | 1.039881e+06 | 9.614689e+05 | ... | 0.000000e+00 | 0.000000e+00 | 4.333662e+05 | 8.622547e+06 | 1.203022e+07 | 15542.090769 | True | 35958 | 12192.0 | 80713.729357 |
| 30 | Upper Harbour Local Board Area | MULTIPOLYGON (((1743880.554 5928979.708, 17438... | 174.672015 | -36.760234 | 6.973017e+07 | 2.098467e+07 | 4.506156e+06 | 7.750997e+06 | 3.291655e+06 | 3.947975e+05 | ... | 9.573434e+06 | 0.000000e+00 | 6.433695e+06 | 1.217632e+07 | 8.853710e+06 | 12662.494326 | True | 53571 | 17079.0 | 92354.882984 |
| 31 | Waiheke | MULTIPOLYGON (((1793981.864 5931917.843, 17940... | 175.055216 | -36.799894 | 1.547768e+08 | 0.000000e+00 | 0.000000e+00 | 0.000000e+00 | 0.000000e+00 | 0.000000e+00 | ... | 0.000000e+00 | 1.521384e+08 | 0.000000e+00 | 0.000000e+00 | 0.000000e+00 | 26689.652804 | True | 8310 | 3597.0 | 55066.946309 |
| 32 | Waitakere | MULTIPOLYGON (((1743890.310 5905057.188, 17438... | 174.527874 | -36.931281 | 2.651435e+08 | 3.686649e+06 | 4.685951e+05 | 1.370472e+05 | 2.243794e+05 | 0.000000e+00 | ... | 2.606228e+05 | 0.000000e+00 | 5.718099e+04 | 8.375348e+07 | 1.672914e+08 | 22820.718458 | True | 12477 | 4449.0 | 84773.121192 |
| 33 | Whau | MULTIPOLYGON (((1748168.244 5916597.058, 17481... | 174.684472 | -36.906902 | 2.682152e+07 | 1.602986e+07 | 9.520285e+05 | 6.562767e+06 | 6.248986e+06 | 2.009848e+06 | ... | 0.000000e+00 | 0.000000e+00 | 3.309943e+06 | 0.000000e+00 | 2.842184e+06 | 9497.105591 | True | 72585 | 23928.0 | 65948.529781 |
33 rows × 21 columns
df.drop(columns="geometry").to_csv("output/Local_Area.csv")
mb = gpd.read_file("input/statsnzmeshblock-higher-geographies-2018-generalised-FGDB.zip!meshblock-higher-geographies-2018-generalised.gdb")
mb = mb[(mb.TA2018_V1_00_NAME == "Auckland") & mb.LANDWATER_NAME.isin(["Mainland", "Island"])]
mb = mb.rename(columns={"MB2018_V1_00": "Code"})
mb["Centroid_lon"] = mb.centroid.to_crs(epsg=4326).x
mb["Centroid_lat"] = mb.centroid.to_crs(epsg=4326).y
mb
| Code | SA12018_V1_00 | SA22018_V1_00 | SA22018_V1_00_NAME | UR2018_V1_00 | UR2018_V1_00_NAME | IUR2018_V1_00 | IUR2018_V1_00_NAME | CB2018_V1_00 | CB2018_V1_00_NAME | ... | WARD2018_V1_00 | WARD2018_V1_00_NAME | LANDWATER | LANDWATER_NAME | LAND_AREA_SQ_KM | AREA_SQ_KM | Shape_Length | geometry | Centroid_lon | Centroid_lat | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 543 | 0137200 | 7001157 | 110400 | Cape Rodney | 1098 | Other rural Auckland | 22 | Rural other | 07601 | Rodney Local Board Area | ... | 07601 | Rodney Ward | 11 | Island | 0.097567 | 0.097567 | 1356.404986 | MULTIPOLYGON (((1761489.925 5985285.678, 17614... | 174.797456 | -36.265557 |
| 544 | 0170700 | 7001317 | 114300 | Gulf Islands | 1098 | Other rural Auckland | 22 | Rural other | 07602 | Hibiscus and Bays Local Board Area | ... | 07602 | Albany Ward | 11 | Island | 0.005545 | 0.005545 | 321.185380 | MULTIPOLYGON (((1753753.763 5954256.206, 17537... | 174.718386 | -36.545155 |
| 545 | 0438700 | 7001130 | 111800 | Barrier Islands | 1098 | Other rural Auckland | 22 | Rural other | 07608 | Great Barrier Local Board Area | ... | 07605 | Waitemata and Gulf Ward | 11 | Island | 70.725594 | 70.725594 | 100201.649851 | MULTIPOLYGON (((1815954.125 6007877.906, 18159... | 175.382085 | -36.108549 |
| 546 | 0439306 | 7001135 | 111800 | Barrier Islands | 1109 | Tryphena | 21 | Rural settlement | 07608 | Great Barrier Local Board Area | ... | 07605 | Waitemata and Gulf Ward | 11 | Island | 0.078233 | 0.078233 | 1494.370224 | MULTIPOLYGON (((1823569.424 5979782.392, 18235... | 175.489541 | -36.302663 |
| 547 | 0439307 | 7001135 | 111800 | Barrier Islands | 1109 | Tryphena | 21 | Rural settlement | 07608 | Great Barrier Local Board Area | ... | 07605 | Waitemata and Gulf Ward | 11 | Island | 0.475337 | 0.475337 | 3408.621475 | MULTIPOLYGON (((1824128.365 5979223.530, 18242... | 175.494343 | -36.309033 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 53299 | 4011884 | 7001155 | 110400 | Cape Rodney | 1073 | Whangateau | 21 | Rural settlement | 07601 | Rodney Local Board Area | ... | 07601 | Rodney Ward | 12 | Mainland | 1.161280 | 1.161280 | 7266.530123 | MULTIPOLYGON (((1759217.632 5980700.819, 17596... | 174.769327 | -36.310203 |
| 53301 | 4011883 | 7001154 | 110400 | Cape Rodney | 1098 | Other rural Auckland | 22 | Rural other | 07601 | Rodney Local Board Area | ... | 07601 | Rodney Ward | 12 | Mainland | 7.680889 | 7.680889 | 17076.458542 | MULTIPOLYGON (((1756196.931 5982761.262, 17562... | 174.752390 | -36.303185 |
| 53333 | 4011925 | 7009199 | 156900 | Baverstock | 1108 | Auckland | 11 | Major urban area | 07616 | Howick Local Board Area | ... | 07610 | Howick Ward | 12 | Mainland | 0.242189 | 0.242189 | 2352.598074 | MULTIPOLYGON (((1770681.901 5909197.803, 17711... | 174.920613 | -36.949951 |
| 53359 | 4011971 | 7009313 | 158600 | Ormiston East | 1108 | Auckland | 11 | Major urban area | 07616 | Howick Local Board Area | ... | 07610 | Howick Ward | 12 | Mainland | 0.487587 | 0.487587 | 3692.342983 | MULTIPOLYGON (((1771223.196 5906518.069, 17715... | 174.930444 | -36.975090 |
| 53360 | 4011972 | 7009335 | 158900 | Tuscany Heights | 1108 | Auckland | 11 | Major urban area | 07616 | Howick Local Board Area | ... | 07610 | Howick Ward | 12 | Mainland | 0.178460 | 0.178460 | 1798.203501 | MULTIPOLYGON (((1772642.183 5906246.394, 17726... | 174.940243 | -36.976625 |
13441 rows × 30 columns
mb = mb.set_index("Code")
skytower = Point(1757109.809, 5920500.841)
mb["Hdist_skytower"] = mb.centroid.distance(skytower) # meters
mb
| SA12018_V1_00 | SA22018_V1_00 | SA22018_V1_00_NAME | UR2018_V1_00 | UR2018_V1_00_NAME | IUR2018_V1_00 | IUR2018_V1_00_NAME | CB2018_V1_00 | CB2018_V1_00_NAME | CON2018_V1_00 | ... | WARD2018_V1_00_NAME | LANDWATER | LANDWATER_NAME | LAND_AREA_SQ_KM | AREA_SQ_KM | Shape_Length | geometry | Centroid_lon | Centroid_lat | Hdist_skytower | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Code | |||||||||||||||||||||
| 0137200 | 7001157 | 110400 | Cape Rodney | 1098 | Other rural Auckland | 22 | Rural other | 07601 | Rodney Local Board Area | 0299 | ... | Rodney Ward | 11 | Island | 0.097567 | 0.097567 | 1356.404986 | MULTIPOLYGON (((1761489.925 5985285.678, 17614... | 174.797456 | -36.265557 | 64744.057618 |
| 0170700 | 7001317 | 114300 | Gulf Islands | 1098 | Other rural Auckland | 22 | Rural other | 07602 | Hibiscus and Bays Local Board Area | 0299 | ... | Albany Ward | 11 | Island | 0.005545 | 0.005545 | 321.185380 | MULTIPOLYGON (((1753753.763 5954256.206, 17537... | 174.718386 | -36.545155 | 33870.034882 |
| 0438700 | 7001130 | 111800 | Barrier Islands | 1098 | Other rural Auckland | 22 | Rural other | 07608 | Great Barrier Local Board Area | 0299 | ... | Waitemata and Gulf Ward | 11 | Island | 70.725594 | 70.725594 | 100201.649851 | MULTIPOLYGON (((1815954.125 6007877.906, 18159... | 175.382085 | -36.108549 | 99126.852773 |
| 0439306 | 7001135 | 111800 | Barrier Islands | 1109 | Tryphena | 21 | Rural settlement | 07608 | Great Barrier Local Board Area | 0299 | ... | Waitemata and Gulf Ward | 11 | Island | 0.078233 | 0.078233 | 1494.370224 | MULTIPOLYGON (((1823569.424 5979782.392, 18235... | 175.489541 | -36.302663 | 88914.804842 |
| 0439307 | 7001135 | 111800 | Barrier Islands | 1109 | Tryphena | 21 | Rural settlement | 07608 | Great Barrier Local Board Area | 0299 | ... | Waitemata and Gulf Ward | 11 | Island | 0.475337 | 0.475337 | 3408.621475 | MULTIPOLYGON (((1824128.365 5979223.530, 18242... | 175.494343 | -36.309033 | 88749.923755 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 4011884 | 7001155 | 110400 | Cape Rodney | 1073 | Whangateau | 21 | Rural settlement | 07601 | Rodney Local Board Area | 0299 | ... | Rodney Ward | 12 | Mainland | 1.161280 | 1.161280 | 7266.530123 | MULTIPOLYGON (((1759217.632 5980700.819, 17596... | 174.769327 | -36.310203 | 59716.438951 |
| 4011883 | 7001154 | 110400 | Cape Rodney | 1098 | Other rural Auckland | 22 | Rural other | 07601 | Rodney Local Board Area | 0299 | ... | Rodney Ward | 12 | Mainland | 7.680889 | 7.680889 | 17076.458542 | MULTIPOLYGON (((1756196.931 5982761.262, 17562... | 174.752390 | -36.303185 | 60497.860816 |
| 4011925 | 7009199 | 156900 | Baverstock | 1108 | Auckland | 11 | Major urban area | 07616 | Howick Local Board Area | 0299 | ... | Howick Ward | 12 | Mainland | 0.242189 | 0.242189 | 2352.598074 | MULTIPOLYGON (((1770681.901 5909197.803, 17711... | 174.920613 | -36.949951 | 18068.596395 |
| 4011971 | 7009313 | 158600 | Ormiston East | 1108 | Auckland | 11 | Major urban area | 07616 | Howick Local Board Area | 0299 | ... | Howick Ward | 12 | Mainland | 0.487587 | 0.487587 | 3692.342983 | MULTIPOLYGON (((1771223.196 5906518.069, 17715... | 174.930444 | -36.975090 | 20557.150466 |
| 4011972 | 7009335 | 158900 | Tuscany Heights | 1108 | Auckland | 11 | Major urban area | 07616 | Howick Local Board Area | 0299 | ... | Howick Ward | 12 | Mainland | 0.178460 | 0.178460 | 1798.203501 | MULTIPOLYGON (((1772642.183 5906246.394, 17726... | 174.940243 | -36.976625 | 21315.702128 |
13441 rows × 30 columns
mb[mb.Hdist_skytower < 5000].plot(column="Hdist_skytower", legend=True)
plt.title("Meshblocks <5KM from the Skytower")
Text(0.5, 1.0, 'Meshblocks <5KM from the Skytower')
coastline = gpd.read_file("input/lds-nz-coastlines-and-islands-polygons-topo-150k-FGDB.zip!nz-coastlines-and-islands-polygons-topo-150k.gdb")
%%time
coastline_shape = coastline.dissolve().boundary.iloc[0]
CPU times: user 34.9 s, sys: 18.9 ms, total: 34.9 s Wall time: 34.9 s
%%time
mb["Hdist_coast"] = mb.centroid.distance(coastline_shape)
CPU times: user 3min 5s, sys: 17.1 ms, total: 3min 5s Wall time: 3min 5s
mb["Hdist_coast"].describe() # meters
count 13441.000000 mean 1596.010824 std 1792.477098 min 0.488582 25% 444.163564 50% 1076.033671 75% 2116.959738 max 16027.664069 Name: Hdist_coast, dtype: float64
ax = mb.plot(column="Hdist_coast", legend=True)
gpd.clip(coastline, mb.dissolve().envelope).boundary.plot(ax=ax)
<AxesSubplot:>
OSRM API docs - http://project-osrm.org/docs/v5.24.0/api/#table-service
len(mb)
13441
BASE_URL = "http://osrm.auckland-cer.cloud.edu.au"
def drive(points, to):
points = [f"{point.x},{point.y}" for point in points]
points = ";".join(points)
result = requests.get(f"{BASE_URL}/table/v1/driving/{to};{points}?destinations=0&annotations=duration,distance")
return result.json()
skytower = "174.76218883819053,-36.848429166610735"
%%time
result = drive(points=mb.centroid.to_crs(epsg=4326), to=skytower)
CPU times: user 1.29 s, sys: 251 µs, total: 1.29 s Wall time: 1min 13s
result.keys()
dict_keys(['code', 'sources', 'destinations', 'durations', 'distances'])
mb["Mdist_skytower"] = [r[0] for r in result["distances"][1:]]
mb["Mdist_skytower"].describe() # Units are meters
count 13441.000000 mean 21083.293989 std 15682.240482 min 59.500000 25% 11506.700000 50% 17201.300000 75% 24218.100000 max 136633.200000 Name: Mdist_skytower, dtype: float64
mb.plot(column="Mdist_skytower", legend=True)
<AxesSubplot:>
mb[mb.Mdist_skytower < 5000].plot(column="Mdist_skytower", legend=True)
plt.title("Meshblocks <5KM driving from the Skytower")
Text(0.5, 1.0, 'Meshblocks <5KM driving from the Skytower')
# Some strange gaps there - let's investigate
points = mb[(mb.Hdist_skytower < 1000) & (mb.Mdist_skytower > 5000)]
points = [f"{point.x},{point.y}" for point in points.centroid.to_crs(epsg=4326)]
points
['174.75624867471396,-36.85383185528419', '174.75345762242358,-36.850503483712096', '174.75551843174853,-36.8526249602678', '174.75734772557055,-36.855642008978414']
routes = []
for point in points:
route_result = requests.get(f"{BASE_URL}/route/v1/driving/{skytower};{point}?geometries=geojson")
routes.extend(route_result.json()["routes"])
routes
[{'geometry': {'coordinates': [[174.762021, -36.848507],
[174.762081, -36.848111],
[174.761663, -36.847994],
[174.760279, -36.851092],
[174.756401, -36.856918],
[174.756298, -36.858244],
[174.757545, -36.859416],
[174.761092, -36.859842],
[174.763064, -36.860822],
[174.764566, -36.862759],
[174.766927, -36.864565],
[174.768415, -36.867436],
[174.774277, -36.87211],
[174.774115, -36.872578],
[174.773385, -36.872143],
[174.771509, -36.870028],
[174.768209, -36.867454],
[174.766746, -36.864588],
[174.764693, -36.863068],
[174.76298, -36.861326],
[174.761243, -36.860494],
[174.758161, -36.860144],
[174.756868, -36.859449],
[174.755729, -36.85724],
[174.756435, -36.854967],
[174.756232, -36.853838]],
'type': 'LineString'},
'legs': [{'steps': [],
'distance': 6623.4,
'duration': 471.2,
'summary': '',
'weight': 482.5}],
'distance': 6623.4,
'duration': 471.2,
'weight_name': 'routability',
'weight': 482.5},
{'geometry': {'coordinates': [[174.762021, -36.848507],
[174.761973, -36.848406],
[174.762081, -36.848111],
[174.761663, -36.847994],
[174.760403, -36.850878],
[174.75828, -36.854101],
[174.758023, -36.853973],
[174.756934, -36.853677],
[174.756692, -36.853485],
[174.756478, -36.853123],
[174.756232, -36.85299],
[174.755149, -36.853046],
[174.75503, -36.853003],
[174.754968, -36.852487],
[174.754131, -36.851451],
[174.753566, -36.85047]],
'type': 'LineString'},
'legs': [{'steps': [],
'distance': 1487.9,
'duration': 201.1,
'summary': '',
'weight': 212.4}],
'distance': 1487.9,
'duration': 201.1,
'weight_name': 'routability',
'weight': 212.4},
{'geometry': {'coordinates': [[174.762021, -36.848507],
[174.762081, -36.848111],
[174.760362, -36.847639],
[174.756384, -36.84814],
[174.75646, -36.846119],
[174.751793, -36.844936],
[174.749101, -36.843285],
[174.746338, -36.842436],
[174.744751, -36.840654],
[174.742077, -36.838893],
[174.741435, -36.837888],
[174.74165, -36.836424],
[174.749967, -36.823771],
[174.750195, -36.822768],
[174.749701, -36.81919],
[174.751533, -36.813923],
[174.750898, -36.812294],
[174.753046, -36.811822],
[174.753968, -36.812071],
[174.754297, -36.812877],
[174.753791, -36.813549],
[174.75206, -36.814361],
[174.751215, -36.815353],
[174.750131, -36.819403],
[174.750332, -36.822928],
[174.750102, -36.823818],
[174.741981, -36.836162],
[174.74164, -36.837484],
[174.742147, -36.838632],
[174.744931, -36.840513],
[174.746304, -36.842152],
[174.749683, -36.843296],
[174.752518, -36.845063],
[174.752984, -36.845821],
[174.753572, -36.850065],
[174.755481, -36.852645]],
'type': 'LineString'},
'legs': [{'steps': [],
'distance': 11247.2,
'duration': 783.1,
'summary': '',
'weight': 794.4}],
'distance': 11247.2,
'duration': 783.1,
'weight_name': 'routability',
'weight': 794.4},
{'geometry': {'coordinates': [[174.762021, -36.848507],
[174.761973, -36.848406],
[174.762081, -36.848111],
[174.761663, -36.847994],
[174.760403, -36.850878],
[174.758622, -36.85365],
[174.757737, -36.854764],
[174.757261, -36.855607]],
'type': 'LineString'},
'legs': [{'steps': [],
'distance': 1022.2,
'duration': 135,
'summary': '',
'weight': 146.3}],
'distance': 1022.2,
'duration': 135,
'weight_name': 'routability',
'weight': 146.3}]
route_df = pd.read_json(json.dumps(routes))
route_df["geometry"] = route_df['geometry'].apply(shape)
route_df = gpd.GeoDataFrame(route_df)
route_df
| geometry | legs | distance | duration | weight_name | weight | |
|---|---|---|---|---|---|---|
| 0 | LINESTRING (174.76202 -36.84851, 174.76208 -36... | [{'steps': [], 'distance': 6623.4, 'duration':... | 6623.4 | 471.2 | routability | 482.5 |
| 1 | LINESTRING (174.76202 -36.84851, 174.76197 -36... | [{'steps': [], 'distance': 1487.9, 'duration':... | 1487.9 | 201.1 | routability | 212.4 |
| 2 | LINESTRING (174.76202 -36.84851, 174.76208 -36... | [{'steps': [], 'distance': 11247.2, 'duration'... | 11247.2 | 783.1 | routability | 794.4 |
| 3 | LINESTRING (174.76202 -36.84851, 174.76197 -36... | [{'steps': [], 'distance': 1022.2, 'duration':... | 1022.2 | 135.0 | routability | 146.3 |
ax = mb[mb.Mdist_skytower < 5000].to_crs(epsg=4326).plot(column="Mdist_skytower", legend=True)
route_df[route_df["distance"]>10000].plot(column="distance", legend=True, ax=ax, linewidth=3, cmap="plasma")
<AxesSubplot:>
points = mb[(mb.Hdist_skytower < 1000) & (mb.Mdist_skytower > 5000)]
points = [f"{point.y},{point.x}" for point in points.centroid.to_crs(epsg=4326)]
for point in points:
print(f"https://map.project-osrm.org/?loc={point}&loc=-36.848429166610735,174.76218883819053&srv=0")
https://map.project-osrm.org/?loc=-36.85383185528419,174.75624867471396&loc=-36.848429166610735,174.76218883819053&srv=0 https://map.project-osrm.org/?loc=-36.850503483712096,174.75345762242358&loc=-36.848429166610735,174.76218883819053&srv=0 https://map.project-osrm.org/?loc=-36.8526249602678,174.75551843174853&loc=-36.848429166610735,174.76218883819053&srv=0 https://map.project-osrm.org/?loc=-36.855642008978414,174.75734772557055&loc=-36.848429166610735,174.76218883819053&srv=0
%%time
north = drive(points=mb.centroid.translate(yoff=-100).to_crs(epsg=4326), to=skytower)
west = drive(points=mb.centroid.translate(xoff=-100).to_crs(epsg=4326), to=skytower)
east = drive(points=mb.centroid.translate(xoff=100).to_crs(epsg=4326), to=skytower)
south = drive(points=mb.centroid.translate(yoff=100).to_crs(epsg=4326), to=skytower)
CPU times: user 6.89 s, sys: 41.7 ms, total: 6.93 s Wall time: 4min 52s
mb["Mdist_skytower_center"] = [r[0] for r in result["distances"][1:]]
mb["Mdist_skytower_north"] = [r[0] for r in north["distances"][1:]]
mb["Mdist_skytower_west"] = [r[0] for r in west["distances"][1:]]
mb["Mdist_skytower_east"] = [r[0] for r in east["distances"][1:]]
mb["Mdist_skytower_south"] = [r[0] for r in south["distances"][1:]]
mb[["Mdist_skytower_center",
"Mdist_skytower_north",
"Mdist_skytower_west",
"Mdist_skytower_east",
"Mdist_skytower_south"]].describe()
| Mdist_skytower_center | Mdist_skytower_north | Mdist_skytower_west | Mdist_skytower_east | Mdist_skytower_south | |
|---|---|---|---|---|---|
| count | 13441.000000 | 13441.000000 | 13441.000000 | 13441.000000 | 13441.000000 |
| mean | 21083.293989 | 21104.614404 | 21076.932974 | 21080.920289 | 21044.649126 |
| std | 15682.240482 | 15672.937314 | 15680.363394 | 15684.086849 | 15692.827461 |
| min | 59.500000 | 29.900000 | 110.700000 | 155.200000 | 109.600000 |
| 25% | 11506.700000 | 11517.900000 | 11485.800000 | 11514.800000 | 11467.600000 |
| 50% | 17201.300000 | 17168.500000 | 17202.300000 | 17154.200000 | 17146.600000 |
| 75% | 24218.100000 | 24276.600000 | 24203.500000 | 24235.700000 | 24178.300000 |
| max | 136633.200000 | 136633.200000 | 136633.200000 | 136633.200000 | 136633.200000 |
mb["Mdist_skytower"] = mb[[
"Mdist_skytower_center",
"Mdist_skytower_north",
"Mdist_skytower_west",
"Mdist_skytower_east",
"Mdist_skytower_south"
]].min(axis=1)
(mb["Mdist_skytower_center"] - mb["Mdist_skytower"]).describe()
count 13441.000000 mean 270.629648 std 506.630148 min 0.000000 25% 81.500000 50% 131.900000 75% 261.300000 max 12976.900000 dtype: float64
mb[mb.Mdist_skytower < 5000].plot(column="Mdist_skytower", legend=True)
plt.title("Meshblocks <5KM driving from the Skytower")
Text(0.5, 1.0, 'Meshblocks <5KM driving from the Skytower')
mb["Mtime_skytower_center"] = [r[0] for r in result["durations"][1:]]
mb["Mtime_skytower_north"] = [r[0] for r in north["durations"][1:]]
mb["Mtime_skytower_west"] = [r[0] for r in west["durations"][1:]]
mb["Mtime_skytower_east"] = [r[0] for r in east["durations"][1:]]
mb["Mtime_skytower_south"] = [r[0] for r in south["durations"][1:]]
display(mb[[
"Mtime_skytower_center",
"Mtime_skytower_north",
"Mtime_skytower_west",
"Mtime_skytower_east",
"Mtime_skytower_south"
]].describe())
mb["Mtime_skytower"] = mb[[
"Mtime_skytower_center",
"Mtime_skytower_north",
"Mtime_skytower_west",
"Mtime_skytower_east",
"Mtime_skytower_south"
]].min(axis=1)
mb["Mtime_skytower"].describe() # Units are seconds
| Mtime_skytower_center | Mtime_skytower_north | Mtime_skytower_west | Mtime_skytower_east | Mtime_skytower_south | |
|---|---|---|---|---|---|
| count | 13441.000000 | 13441.000000 | 13441.000000 | 13441.000000 | 13441.000000 |
| mean | 1594.367316 | 1593.642854 | 1592.224001 | 1592.324723 | 1590.448925 |
| std | 3358.619612 | 3358.274390 | 3358.655563 | 3357.112345 | 3358.303358 |
| min | 15.100000 | 7.200000 | 26.700000 | 38.900000 | 37.900000 |
| 25% | 867.000000 | 863.700000 | 865.100000 | 862.400000 | 862.400000 |
| 50% | 1163.400000 | 1160.600000 | 1160.800000 | 1160.200000 | 1159.400000 |
| 75% | 1514.500000 | 1510.900000 | 1510.600000 | 1511.700000 | 1508.300000 |
| max | 70521.800000 | 70521.800000 | 70521.800000 | 70521.800000 | 70521.800000 |
count 13441.000000 mean 1568.247095 std 3357.297496 min 7.200000 25% 841.300000 50% 1138.400000 75% 1485.000000 max 70521.800000 Name: Mtime_skytower, dtype: float64
mb[mb.Mtime_skytower < 60*20].plot(column="Mtime_skytower", legend=True)
plt.title("Meshblocks < 20 minutes drive from the Skytower")
Text(0.5, 1.0, 'Meshblocks < 20 minutes drive from the Skytower')
transit = pd.read_excel("input/Geocoordinates_Direct_Transit_stops_AKL.xlsx")
transit = gpd.GeoDataFrame(transit, geometry=gpd.points_from_xy(transit.Longitude, transit.Latitude), crs="EPSG:4326")
transit
| Location | Latitude | Longitude | geometry | |
|---|---|---|---|---|
| 0 | West Harbour | -36.810902 | 174.645571 | POINT (174.64557 -36.81090) |
| 1 | Hobsonville | -36.787563 | 174.672260 | POINT (174.67226 -36.78756) |
| 2 | Beach Haven | -36.789800 | 174.678317 | POINT (174.67832 -36.78980) |
| 3 | Birkenhead | -36.822961 | 174.733882 | POINT (174.73388 -36.82296) |
| 4 | Northcote | -36.826891 | 174.746329 | POINT (174.74633 -36.82689) |
| 5 | Bayswater | -36.821697 | 174.766296 | POINT (174.76630 -36.82170) |
| 6 | Stanley Bay | -36.828009 | 174.781236 | POINT (174.78124 -36.82801) |
| 7 | Downtown | -36.842912 | 174.766930 | POINT (174.76693 -36.84291) |
| 8 | Devonport | -36.833343 | 174.795631 | POINT (174.79563 -36.83334) |
| 9 | Matiatia | -36.780617 | 174.991460 | POINT (174.99146 -36.78062) |
| 10 | Pine Harbour | -36.889582 | 174.989955 | POINT (174.98996 -36.88958) |
| 11 | Half Moon Bay | -36.879913 | 174.897656 | POINT (174.89766 -36.87991) |
| 12 | Gulf Harbour | -36.624279 | 174.787779 | POINT (174.78778 -36.62428) |
| 13 | Avondale | -36.897302 | 174.699092 | POINT (174.69909 -36.89730) |
| 14 | Baldwin Avenue | -36.877694 | 174.720462 | POINT (174.72046 -36.87769) |
| 15 | Britomart | -36.844177 | 174.767766 | POINT (174.76777 -36.84418) |
| 16 | Ellerslie | -36.898225 | 174.808085 | POINT (174.80809 -36.89822) |
| 17 | Fruitvale Road | -36.910636 | 174.667006 | POINT (174.66701 -36.91064) |
| 18 | Glen Eden | -36.910149 | 174.652854 | POINT (174.65285 -36.91015) |
| 19 | Glen Innes | -36.878788 | 174.854118 | POINT (174.85412 -36.87879) |
| 20 | Grafton | -36.865501 | 174.770008 | POINT (174.77001 -36.86550) |
| 21 | Greenlane | -36.889649 | 174.797444 | POINT (174.79744 -36.88965) |
| 22 | Henderson | -36.880903 | 174.630923 | POINT (174.63092 -36.88090) |
| 23 | Homai | -37.013442 | 174.874655 | POINT (174.87465 -37.01344) |
| 24 | Kingsland | -36.872443 | 174.744520 | POINT (174.74452 -36.87244) |
| 25 | Manukau | -36.993769 | 174.877389 | POINT (174.87739 -36.99377) |
| 26 | Manurewa | -37.023261 | 174.896139 | POINT (174.89614 -37.02326) |
| 27 | Meadowbank | -36.866298 | 174.820777 | POINT (174.82078 -36.86630) |
| 28 | Middlemore | -36.963044 | 174.838964 | POINT (174.83896 -36.96304) |
| 29 | Morningside | -36.874913 | 174.735210 | POINT (174.73521 -36.87491) |
| 30 | Mt Eden | -36.867955 | 174.758970 | POINT (174.75897 -36.86796) |
| 31 | Mt. Albert | -36.884793 | 174.714007 | POINT (174.71401 -36.88479) |
| 32 | New Lynn | -36.909357 | 174.684020 | POINT (174.68402 -36.90936) |
| 33 | Newmarket | -36.869622 | 174.778879 | POINT (174.77888 -36.86962) |
| 34 | Onehunga | -36.925383 | 174.786012 | POINT (174.78601 -36.92538) |
| 35 | Orakei | -36.862352 | 174.809466 | POINT (174.80947 -36.86235) |
| 36 | Otahuhu | -36.946888 | 174.833261 | POINT (174.83326 -36.94689) |
| 37 | Panmure | -36.898136 | 174.849269 | POINT (174.84927 -36.89814) |
| 38 | Papatoetoe | -36.977583 | 174.849376 | POINT (174.84938 -36.97758) |
| 39 | Papakura | -37.064926 | 174.946311 | POINT (174.94631 -37.06493) |
| 40 | Parnell | -36.854689 | 174.777457 | POINT (174.77746 -36.85469) |
| 41 | Penrose | -36.910093 | 174.815687 | POINT (174.81569 -36.91009) |
| 42 | Puhinui | -36.989788 | 174.856028 | POINT (174.85603 -36.98979) |
| 43 | Pukekohe | -37.203246 | 174.910155 | POINT (174.91016 -37.20325) |
| 44 | Ranui | -36.867874 | 174.603303 | POINT (174.60330 -36.86787) |
| 45 | Remuera | -36.881313 | 174.785404 | POINT (174.78540 -36.88131) |
| 46 | Sturges Road | -36.873433 | 174.620874 | POINT (174.62087 -36.87343) |
| 47 | Sunnyvale | -36.896779 | 174.631982 | POINT (174.63198 -36.89678) |
| 48 | Swanson | -36.866032 | 174.576258 | POINT (174.57626 -36.86603) |
| 49 | Sylvia Park | -36.914619 | 174.842589 | POINT (174.84259 -36.91462) |
| 50 | Takanini | -37.041127 | 174.919446 | POINT (174.91945 -37.04113) |
| 51 | Te Mahia | -37.031118 | 174.906152 | POINT (174.90615 -37.03112) |
| 52 | Te Papapa | -36.920079 | 174.801432 | POINT (174.80143 -36.92008) |
| 53 | Hibiscus Coast Bus Station | -36.624647 | 174.667158 | POINT (174.66716 -36.62465) |
| 54 | Constellation Bus Station | -36.752014 | 174.728166 | POINT (174.72817 -36.75201) |
| 55 | Sunnynook Bus Station | -36.761481 | 174.737999 | POINT (174.73800 -36.76148) |
| 56 | Smales Farm Bus Station | -36.784523 | 174.751077 | POINT (174.75108 -36.78452) |
| 57 | Akoranga Bus Station | -36.797232 | 174.760938 | POINT (174.76094 -36.79723) |
len(mb)
13441
%%time
def drive(from_points, to_points):
# Drive distance and time these meshblock centroids to all 58 rapid transit stops
from_points_s = ";".join([f"{point.x},{point.y}" for point in from_points])
from_indices = ";".join([str(i) for i in range(len(from_points))])
to_points_s = ";".join([f"{point.x},{point.y}" for point in to_points])
to_indices = ";".join([str(i + len(from_points)) for i in range(len(to_points))])
result = requests.get(f"{BASE_URL}/table/v1/driving/{from_points_s};{to_points_s}?sources={from_indices}&destinations={to_indices}&annotations=duration,distance")
return result.json()
result = drive(from_points=mb.centroid.to_crs(epsg=4326), to_points=transit.geometry.to_crs(epsg=4326))
CPU times: user 1.77 s, sys: 39.2 ms, total: 1.81 s Wall time: 30.7 s
min_indices = np.argmin(result["distances"], axis=1)
mb["Mdist_rapid_transit"] = [result["distances"][i][min_indices[i]] for i in range(len(min_indices))]
mb["Mtime_rapid_transit"] = [result["durations"][i][min_indices[i]] for i in range(len(min_indices))]
mb["Mdist_rapid_transit_name"] = list(transit.Location[min_indices])
mb.Mdist_rapid_transit_name.value_counts()
Constellation Bus Station 1006 Hibiscus Coast Bus Station 862 Half Moon Bay 718 Pukekohe 515 Manukau 503 Papakura 499 Glen Innes 424 Onehunga 420 Manurewa 378 Kingsland 371 Mt. Albert 346 West Harbour 332 Glen Eden 326 Middlemore 314 Avondale 304 Henderson 288 Papatoetoe 286 Fruitvale Road 248 Britomart 246 Otahuhu 246 Morningside 227 Beach Haven 223 Akoranga Bus Station 223 Sunnyvale 219 Sylvia Park 218 Mt Eden 216 Smales Farm Bus Station 214 Ranui 208 New Lynn 202 Gulf Harbour 193 Newmarket 185 Orakei 182 Sturges Road 177 Homai 154 Birkenhead 147 Matiatia 144 Takanini 128 Ellerslie 125 Panmure 122 Swanson 121 Meadowbank 118 Devonport 115 Grafton 101 Te Papapa 96 Bayswater 95 Remuera 91 Pine Harbour 90 Parnell 90 Puhinui 83 Penrose 81 Baldwin Avenue 68 Northcote 68 Hobsonville 56 Stanley Bay 15 Te Mahia 11 Sunnynook Bus Station 2 Downtown 1 Name: Mdist_rapid_transit_name, dtype: int64
ax = mb.to_crs(transit.crs).plot(column="Mdist_rapid_transit", legend=True)
transit.plot(ax=ax, color="red")
<AxesSubplot:>
ax = mb[mb.Mdist_rapid_transit < 5000].to_crs(transit.crs).plot(column="Mdist_rapid_transit", legend=True)
transit.plot(ax=ax, color="red")
<AxesSubplot:>
onramps = gpd.read_file("input/kx-nz-state-highway-on-ramps-off-ramps-SHP.zip")
onramps = gpd.clip(onramps.to_crs(mb.crs), mb)
onramps
| ID | ROADID | SH | INTERCHANG | DISPLACEME | RAMPNO | LENGTH | DESCRIPTIO | geometry | |
|---|---|---|---|---|---|---|---|---|---|
| 1 | 288.0 | 288 | 020 | 0019 | 0.00 | 1 | 37.0 | 020-0010/03.65-D-ON | LINESTRING (1758742.502 5911581.839, 1758734.3... |
| 2 | 251.0 | 251 | 016 | 0027 | 0.00 | 1 | 155.0 | 016-0007/04.01-I-OFF | LINESTRING (1747446.038 5919555.363, 1747389.5... |
| 5 | 291.0 | 291 | 020 | 0022 | 0.00 | 1 | 203.0 | 020-0010/02.28-I-ON | LINESTRING (1759420.418 5910233.008, 1759414.2... |
| 34 | 1765.0 | 1765 | 020 | 0002 | 0.00 | 1 | 414.0 | 020-0010/04.77-I-OFF | LINESTRING (1757782.109 5912150.461, 1757752.3... |
| 35 | 2623.0 | 2623 | 01N | 2036 | 0.00 | 1 | 337.0 | 01N-2398/12.36-X410-R3-OFF | LINESTRING (1753039.078 5934634.671, 1753012.1... |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 437 | 6297867.0 | 3101 | 018 | 0005 | None | 1 | 313.0 | 018-0005-R1 | LINESTRING (1750624.813 5929454.552, 1750604.6... |
| 438 | 6297864.0 | 3106 | 018 | 0002 | None | 4 | 263.0 | 018-0002-R4 | LINESTRING (1752046.399 5930728.423, 1752081.5... |
| 439 | 6297838.0 | 3086 | 016 | 0100 | None | 2 | 843.0 | 016-0003-R2 | LINESTRING (1756829.783 5919210.950, 1756778.7... |
| 440 | 6297836.0 | 3088 | 016 | 0002 | None | 4 | 437.0 | 016-0002-R4 TEMP ON | LINESTRING (1757567.004 5919997.142, 1757604.0... |
| 446 | 6298167.0 | 3253 | 016 | 0019 | None | 3 | 68.0 | 016-0019-R3 | LINESTRING (1743974.769 5923934.492, 1743962.0... |
256 rows × 9 columns
# Filter out offramps
onramps = onramps[~onramps.DESCRIPTIO.str.contains("OFF")]
onramps
| ID | ROADID | SH | INTERCHANG | DISPLACEME | RAMPNO | LENGTH | DESCRIPTIO | geometry | |
|---|---|---|---|---|---|---|---|---|---|
| 1 | 288.0 | 288 | 020 | 0019 | 0.00 | 1 | 37.0 | 020-0010/03.65-D-ON | LINESTRING (1758742.502 5911581.839, 1758734.3... |
| 5 | 291.0 | 291 | 020 | 0022 | 0.00 | 1 | 203.0 | 020-0010/02.28-I-ON | LINESTRING (1759420.418 5910233.008, 1759414.2... |
| 37 | 2622.0 | 2622 | 01N | 2035 | 0.00 | 1 | 98.0 | 01N-2398/13.96-X412-R4-ON | LINESTRING (1753301.310 5933060.772, 1753303.1... |
| 38 | 258.0 | 258 | 016 | 0034 | 0.00 | 1 | 240.0 | 016-0007/06.43-I-ON | LINESTRING (1745470.285 5920590.583, 1745458.7... |
| 39 | 1981.0 | 1981 | 01N | 2018 | 0.00 | 1 | 352.0 | 01N-2398/15.45-X414-R4-ON | LINESTRING (1754036.875 5931380.043, 1753985.1... |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 437 | 6297867.0 | 3101 | 018 | 0005 | None | 1 | 313.0 | 018-0005-R1 | LINESTRING (1750624.813 5929454.552, 1750604.6... |
| 438 | 6297864.0 | 3106 | 018 | 0002 | None | 4 | 263.0 | 018-0002-R4 | LINESTRING (1752046.399 5930728.423, 1752081.5... |
| 439 | 6297838.0 | 3086 | 016 | 0100 | None | 2 | 843.0 | 016-0003-R2 | LINESTRING (1756829.783 5919210.950, 1756778.7... |
| 440 | 6297836.0 | 3088 | 016 | 0002 | None | 4 | 437.0 | 016-0002-R4 TEMP ON | LINESTRING (1757567.004 5919997.142, 1757604.0... |
| 446 | 6298167.0 | 3253 | 016 | 0019 | None | 3 | 68.0 | 016-0019-R3 | LINESTRING (1743974.769 5923934.492, 1743962.0... |
173 rows × 9 columns
onramps["has_on"] = onramps.DESCRIPTIO.str.contains("ON").astype(str)
onramps.plot(column="has_on", legend=True)
<AxesSubplot:>
%%time
result = drive(from_points=mb.centroid.to_crs(epsg=4326), to_points=onramps.centroid.to_crs(epsg=4326))
CPU times: user 2.34 s, sys: 82.1 ms, total: 2.42 s Wall time: 1min
len(result["distances"]), len(result["distances"][0])
(13441, 173)
min_indices = np.argmin(result["distances"], axis=1)
mb["Mdist_onramp"] = [result["distances"][i][min_indices[i]] for i in range(len(min_indices))]
mb["Mtime_onramp"] = [result["durations"][i][min_indices[i]] for i in range(len(min_indices))]
mb["Mdist_onramp_name"] = list(onramps.DESCRIPTIO.iloc[min_indices])
mb.Mdist_onramp_name.value_counts()
01N-2431/07.16-X438-R4-ON 719
016-0007/04.39-I-ON 663
01N-0434-R6 575
01N-2431/13.34-X444-R2-ON 560
016-0007/06.46-I-ON1 487
...
016-0019-R2 1
01N-2398/15.45-X414-R4-ON 1
018-0008-R1 1
01N-0389-R4 1
016-0000/08.53-X8-R4-ON 1
Name: Mdist_onramp_name, Length: 115, dtype: int64
ax = mb.plot(column="Mdist_onramp", legend=True)
onramps.plot(ax=ax, color="red")
<AxesSubplot:>
ax = mb[mb.Mdist_onramp < 5000].plot(column="Mdist_onramp", legend=True)
onramps.plot(ax=ax, color="red")
<AxesSubplot:>
coastline = gpd.read_file("input/lds-nz-coastlines-and-islands-polygons-topo-150k-FGDB.zip!nz-coastlines-and-islands-polygons-topo-150k.gdb")
%%time
coastline = gpd.clip(coastline, mb.dissolve().envelope)
CPU times: user 17.3 s, sys: 9.96 ms, total: 17.3 s Wall time: 17.3 s
%%time
coastline = coastline.dissolve().boundary.buffer(100)
coastline.plot()
CPU times: user 1min 18s, sys: 140 ms, total: 1min 18s Wall time: 1min 18s
<AxesSubplot:>
coastline = coastline.iloc[0]
%%time
mb["Coast_indicator"] = mb.geometry.progress_apply(lambda poly: coastline.intersects(poly))
CPU times: user 16min 54s, sys: 1.29 s, total: 16min 55s Wall time: 16min 52s
mb.plot(column="Coast_indicator", legend=True)
<AxesSubplot:>
# Note that this is meshblock 2020, which is not ideal, but there doesn't seem to be any changes in the Auckland region between 2018 and 2020, so it's probably fine
pop = pd.read_csv("input/2018-census-electoral-population-meshblock-2020-data.csv", index_col=0)
pop
| General_Electoral_Population | Maori_Electoral_Population | GED2020_V1_00 | GED2020_V1_00_NAME | GED2020_V1_00_NAME_ASCII | MED2020_V1_00 | MED2020_V1_00_NAME | MED2020_V1_00_NAME_ASCII | LAND_AREA_SQ_KM | AREA_SQ_KM | |
|---|---|---|---|---|---|---|---|---|---|---|
| MB2020_V2_00 | ||||||||||
| 100 | -999 | 9 | 32 | Northland | Northland | 5 | Te Tai Tokerau | Te Tai Tokerau | 157.497825 | 157.497825 |
| 200 | 9 | 66 | 32 | Northland | Northland | 5 | Te Tai Tokerau | Te Tai Tokerau | 120.503770 | 120.503770 |
| 300 | 12 | 48 | 32 | Northland | Northland | 5 | Te Tai Tokerau | Te Tai Tokerau | 7.481859 | 7.481859 |
| 400 | 6 | 9 | 32 | Northland | Northland | 5 | Te Tai Tokerau | Te Tai Tokerau | 83.342952 | 83.342952 |
| 501 | -999 | -999 | 32 | Northland | Northland | 5 | Te Tai Tokerau | Te Tai Tokerau | 0.000000 | 63.825713 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 4012026 | 24 | -999 | 60 | Wellington Central | Wellington Central | 6 | Te Tai Tonga | Te Tai Tonga | 8.301605 | 8.301605 |
| 4012027 | 108 | 6 | 56 | Waikato | Waikato | 1 | Hauraki-Waikato | Hauraki-Waikato | 5.064385 | 5.064385 |
| 4012028 | 57 | -999 | 51 | Taupō | Taupo | 1 | Hauraki-Waikato | Hauraki-Waikato | 2.612804 | 2.612804 |
| 4012029 | -999 | -999 | 47 | Taieri | Taieri | 6 | Te Tai Tonga | Te Tai Tonga | 0.000000 | 143.908480 |
| 4012030 | -999 | -999 | 8 | Dunedin | Dunedin | 6 | Te Tai Tonga | Te Tai Tonga | 0.000000 | 1016.346178 |
53582 rows × 10 columns
pop.General_Electoral_Population = pop.General_Electoral_Population.replace(-999, np.nan)
mb["Census2018_population"] = list(pop.General_Electoral_Population[mb.index.astype(int)])
display(mb["Census2018_population"].describe())
mb.plot(column="Census2018_population", legend=True)
count 12887.000000 mean 113.744859 std 63.156071 min 6.000000 25% 72.000000 50% 108.000000 75% 147.000000 max 1029.000000 Name: Census2018_population, dtype: float64
<AxesSubplot:>
# Can't do this one due to missing data
mb
| SA12018_V1_00 | SA22018_V1_00 | SA22018_V1_00_NAME | UR2018_V1_00 | UR2018_V1_00_NAME | IUR2018_V1_00 | IUR2018_V1_00_NAME | CB2018_V1_00 | CB2018_V1_00_NAME | CON2018_V1_00 | ... | Mtime_skytower_south | Mtime_skytower | Mdist_rapid_transit | Mtime_rapid_transit | Mdist_rapid_transit_name | Mdist_onramp | Mtime_onramp | Mdist_onramp_name | Coast_indicator | Census2018_population | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Code | |||||||||||||||||||||
| 0137200 | 7001157 | 110400 | Cape Rodney | 1098 | Other rural Auckland | 22 | Rural other | 07601 | Rodney Local Board Area | 0299 | ... | 4741.4 | 4699.9 | 54616.0 | 3348.4 | Hibiscus Coast Bus Station | 43011.9 | 2735.9 | 01N-0381-R2 | True | NaN |
| 0170700 | 7001317 | 114300 | Gulf Islands | 1098 | Other rural Auckland | 22 | Rural other | 07602 | Hibiscus and Bays Local Board Area | 0299 | ... | 2432.5 | 2432.5 | 14338.0 | 1039.5 | Hibiscus Coast Bus Station | 6005.1 | 428.9 | 01N-0381-R2 | True | NaN |
| 0438700 | 7001130 | 111800 | Barrier Islands | 1098 | Other rural Auckland | 22 | Rural other | 07608 | Great Barrier Local Board Area | 0299 | ... | 70462.4 | 70462.4 | 137691.6 | 70755.9 | Britomart | 136507.4 | 70483.9 | 01N-2414/13.17-X427-R4-ON | True | 51.0 |
| 0439306 | 7001135 | 111800 | Barrier Islands | 1109 | Tryphena | 21 | Rural settlement | 07608 | Great Barrier Local Board Area | 0299 | ... | 65818.6 | 65762.6 | 97961.5 | 66062.6 | Britomart | 96777.3 | 65790.6 | 01N-2414/13.17-X427-R4-ON | True | 18.0 |
| 0439307 | 7001135 | 111800 | Barrier Islands | 1109 | Tryphena | 21 | Rural settlement | 07608 | Great Barrier Local Board Area | 0299 | ... | 65724.3 | 65562.3 | 97683.6 | 66021.3 | Britomart | 96499.5 | 65749.3 | 01N-2414/13.17-X427-R4-ON | True | 48.0 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 4011884 | 7001155 | 110400 | Cape Rodney | 1073 | Whangateau | 21 | Rural settlement | 07601 | Rodney Local Board Area | 0299 | ... | 3997.2 | 3987.2 | 46452.6 | 2600.8 | Hibiscus Coast Bus Station | 34848.5 | 1988.3 | 01N-0381-R2 | True | 108.0 |
| 4011883 | 7001154 | 110400 | Cape Rodney | 1098 | Other rural Auckland | 22 | Rural other | 07601 | Rodney Local Board Area | 0299 | ... | 4143.4 | 3983.3 | 46892.7 | 2745.1 | Hibiscus Coast Bus Station | 35288.6 | 2132.6 | 01N-0381-R2 | True | 102.0 |
| 4011925 | 7009199 | 156900 | Baverstock | 1108 | Auckland | 11 | Major urban area | 07616 | Howick Local Board Area | 0299 | ... | 1498.7 | 1486.2 | 7873.1 | 618.7 | Manukau | 6290.5 | 524.9 | 01N-2431/13.34-X444-R2-ON | False | 132.0 |
| 4011971 | 7009313 | 158600 | Ormiston East | 1108 | Auckland | 11 | Major urban area | 07616 | Howick Local Board Area | 0299 | ... | 1580.3 | 1562.8 | 6404.7 | 573.7 | Manukau | 5538.2 | 479.1 | 01N-0448-R2 | False | 81.0 |
| 4011972 | 7009335 | 158900 | Tuscany Heights | 1108 | Auckland | 11 | Major urban area | 07616 | Howick Local Board Area | 0299 | ... | 1744.1 | 1739.0 | 10376.0 | 866.0 | Manukau | 8993.9 | 772.6 | 01N-2431/13.34-X444-R2-ON | False | 66.0 |
13441 rows × 51 columns
mb.drop(columns="geometry").to_csv("output/2018_Meshblocks.csv")
sa2 = gpd.read_file("input/statsnzstatistical-area-2-higher-geographies-2018-clipped-generalis-FGDB.zip!statistical-area-2-higher-geographies-2018-clipped-generalis.gdb")
sa2 = sa2[sa2.TA2018_V1_00_NAME == "Auckland"]
sa2
| SA22018_V1_00 | SA22018_V1_00_NAME | REGC2018_V1_00 | REGC2018_V1_00_NAME | TA2018_V1_00 | TA2018_V1_00_NAME | LAND_AREA_SQ_KM | AREA_SQ_KM | Shape_Length | geometry | |
|---|---|---|---|---|---|---|---|---|---|---|
| 33 | 110900 | Dome Valley-Matakana | 02 | Auckland Region | 076 | Auckland | 84.749966 | 84.749966 | 70031.926112 | MULTIPOLYGON (((1749212.212 5981874.734, 17492... |
| 34 | 111100 | Warkworth West | 02 | Auckland Region | 076 | Auckland | 11.359752 | 11.359752 | 22629.664017 | MULTIPOLYGON (((1748589.106 5974509.939, 17485... |
| 35 | 111200 | Puhoi Valley | 02 | Auckland Region | 076 | Auckland | 236.529965 | 236.529965 | 142134.544465 | MULTIPOLYGON (((1740248.959 5974374.891, 17402... |
| 36 | 112600 | Waikoukou Valley | 02 | Auckland Region | 076 | Auckland | 52.854583 | 52.854583 | 42565.729322 | MULTIPOLYGON (((1735693.112 5942203.564, 17359... |
| 37 | 112800 | Hatfields Beach | 02 | Auckland Region | 076 | Auckland | 0.923255 | 0.923255 | 4314.440531 | MULTIPOLYGON (((1750981.811 5952164.661, 17512... |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 1060 | 165600 | Cape Hill | 02 | Auckland Region | 076 | Auckland | 0.717375 | 0.717375 | 3569.522355 | MULTIPOLYGON (((1770139.188 5882466.221, 17701... |
| 1061 | 165800 | Rooseville Park | 02 | Auckland Region | 076 | Auckland | 1.345809 | 1.345809 | 6010.945714 | MULTIPOLYGON (((1770089.518 5882192.248, 17701... |
| 1062 | 166000 | Pukekohe Central | 02 | Auckland Region | 076 | Auckland | 2.939129 | 2.939129 | 9211.796424 | MULTIPOLYGON (((1769469.158 5881504.407, 17694... |
| 1063 | 166300 | Bombay Hills | 02 | Auckland Region | 076 | Auckland | 31.682020 | 31.682020 | 37482.345738 | MULTIPOLYGON (((1777896.083 5888522.045, 17779... |
| 1064 | 166400 | Ararimu | 02 | Auckland Region | 076 | Auckland | 91.308241 | 91.308241 | 57552.107386 | MULTIPOLYGON (((1786544.313 5891035.131, 17865... |
556 rows × 10 columns
sa2 = sa2.rename(columns={"SA22018_V1_00_NAME": "Name", "SA22018_V1_00": "Code"})
sa2["Centroid_lon"] = sa2.centroid.to_crs(epsg=4326).x
sa2["Centroid_lat"] = sa2.centroid.to_crs(epsg=4326).y
sa2
| Code | Name | REGC2018_V1_00 | REGC2018_V1_00_NAME | TA2018_V1_00 | TA2018_V1_00_NAME | LAND_AREA_SQ_KM | AREA_SQ_KM | Shape_Length | geometry | Centroid_lon | Centroid_lat | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 33 | 110900 | Dome Valley-Matakana | 02 | Auckland Region | 076 | Auckland | 84.749966 | 84.749966 | 70031.926112 | MULTIPOLYGON (((1749212.212 5981874.734, 17492... | 174.654164 | -36.346613 |
| 34 | 111100 | Warkworth West | 02 | Auckland Region | 076 | Auckland | 11.359752 | 11.359752 | 22629.664017 | MULTIPOLYGON (((1748589.106 5974509.939, 17485... | 174.651601 | -36.388932 |
| 35 | 111200 | Puhoi Valley | 02 | Auckland Region | 076 | Auckland | 236.529965 | 236.529965 | 142134.544465 | MULTIPOLYGON (((1740248.959 5974374.891, 17402... | 174.624752 | -36.468509 |
| 36 | 112600 | Waikoukou Valley | 02 | Auckland Region | 076 | Auckland | 52.854583 | 52.854583 | 42565.729322 | MULTIPOLYGON (((1735693.112 5942203.564, 17359... | 174.515741 | -36.709892 |
| 37 | 112800 | Hatfields Beach | 02 | Auckland Region | 076 | Auckland | 0.923255 | 0.923255 | 4314.440531 | MULTIPOLYGON (((1750981.811 5952164.661, 17512... | 174.689356 | -36.570019 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 1060 | 165600 | Cape Hill | 02 | Auckland Region | 076 | Auckland | 0.717375 | 0.717375 | 3569.522355 | MULTIPOLYGON (((1770139.188 5882466.221, 17701... | 174.910713 | -37.187970 |
| 1061 | 165800 | Rooseville Park | 02 | Auckland Region | 076 | Auckland | 1.345809 | 1.345809 | 6010.945714 | MULTIPOLYGON (((1770089.518 5882192.248, 17701... | 174.913076 | -37.196977 |
| 1062 | 166000 | Pukekohe Central | 02 | Auckland Region | 076 | Auckland | 2.939129 | 2.939129 | 9211.796424 | MULTIPOLYGON (((1769469.158 5881504.407, 17694... | 174.913325 | -37.208920 |
| 1063 | 166300 | Bombay Hills | 02 | Auckland Region | 076 | Auckland | 31.682020 | 31.682020 | 37482.345738 | MULTIPOLYGON (((1777896.083 5888522.045, 17779... | 175.002256 | -37.175368 |
| 1064 | 166400 | Ararimu | 02 | Auckland Region | 076 | Auckland | 91.308241 | 91.308241 | 57552.107386 | MULTIPOLYGON (((1786544.313 5891035.131, 17865... | 175.082470 | -37.146862 |
556 rows × 12 columns
sa2["Centroid_local_area"] = sa2.representative_point().apply(lambda point: df.Name[df.contains(point)].iloc[0] if any(df.contains(point)) else np.nan)
sa2["Centroid_local_area"].value_counts()
Henderson-Massey 37 Auckland Central 32 Beach Haven-Birkenhead-Northcote 28 Manurewa 28 Orakei 27 Mangere-Otahuhu 26 Whau 26 Botany 23 Upper Harbour Local Board Area 22 Papakura 22 Hibiscus Coast 22 Puketapapa 20 Devonport-Takapuna 19 Franklin-Pukekohe 17 Mt Albert 16 Tamaki 16 Papatoetoe 16 Mt Eden 15 East Coast Bays 15 Pakuranga 14 Rodney-Helensville 13 Titirangi 12 Howick 12 Maungakiekie 11 Franklin-Beachlands-Hunua 11 Otara 10 Rodney-Warkworth 9 Waitakere 9 Franklin-Waiuku 7 Waiheke 7 Rodney-Kumeu-Riverhead 5 Rodney-Dairy Flat 4 Rodney-Wellsford 4 Name: Centroid_local_area, dtype: int64
sa2[pd.isna(sa2.Centroid_local_area)]
| Code | Name | REGC2018_V1_00 | REGC2018_V1_00_NAME | TA2018_V1_00 | TA2018_V1_00_NAME | LAND_AREA_SQ_KM | AREA_SQ_KM | Shape_Length | geometry | Centroid_lon | Centroid_lat | Centroid_local_area | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 629 | 111800 | Barrier Islands | 02 | Auckland Region | 076 | Auckland | 320.414101 | 320.414101 | 371260.805871 | MULTIPOLYGON (((1815954.125 6007877.906, 18159... | 175.387745 | -36.19825 | NaN |
skytower = Point(1757109.809, 5920500.841) # in NZGD2000 projection
sa2["Hdist_skytower"] = sa2.centroid.distance(skytower)
sa2.plot(column="Hdist_skytower", legend=True)
<AxesSubplot:>
mway = gpd.read_file("input/lds-nz-road-centrelines-topo-150k-FGDB.zip!nz-road-centrelines-topo-150k.gdb")
mway = mway[~pd.isna(mway.hway_num)]
mway
| t50_fid | name_ascii | macronated | name | hway_num | rna_sufi | lane_count | way_count | status | surface | geometry | |
|---|---|---|---|---|---|---|---|---|---|---|---|
| 46 | 100120164 | OHURA ROAD | N | OHURA ROAD | 43 | 3063799 | 2 | None | None | sealed | MULTILINESTRING ((1726407.913 5648240.420, 172... |
| 56 | 100047464 | STATE HIGHWAY 27 | N | STATE HIGHWAY 27 | 27 | 3031731 | 2 | None | None | sealed | MULTILINESTRING ((1844110.982 5802925.948, 184... |
| 62 | 100047472 | STATE HIGHWAY 29A | N | STATE HIGHWAY 29A | 29A | 3025693 | 2 | None | None | sealed | MULTILINESTRING ((1876759.340 5818523.003, 187... |
| 68 | 100047481 | STATE HIGHWAY 2 | N | STATE HIGHWAY 2 | 2 | 3045160 | 2 | None | None | sealed | MULTILINESTRING ((1875605.931 5823398.039, 187... |
| 77 | 100089761 | WEST COAST ROAD | N | WEST COAST ROAD | 73 | 3072232 | 2 | None | None | sealed | MULTILINESTRING ((1490507.178 5234796.972, 149... |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 138626 | 100121061 | STATE HIGHWAY 5 | N | STATE HIGHWAY 5 | 5 | 3068040 | 2 | None | None | sealed | MULTILINESTRING ((1884714.508 5730000.000, 188... |
| 138635 | 100121070 | CHRISTCHURCH SOUTHERN MOTORWAY | N | CHRISTCHURCH SOUTHERN MOTORWAY | 76 | 3061191 | 4 | None | None | sealed | MULTILINESTRING ((1555118.439 5176036.767, 155... |
| 138636 | 100121071 | WESTERN BELFAST BYPASS | N | WESTERN BELFAST BYPASS | 1 | 3079158 | 4 | None | None | sealed | MULTILINESTRING ((1567931.555 5188634.481, 156... |
| 138637 | 100121072 | CHRISTCHURCH NORTHERN MOTORWAY | N | CHRISTCHURCH NORTHERN MOTORWAY | 74 | 0 | 4 | None | None | sealed | MULTILINESTRING ((1571390.312 5190959.558, 157... |
| 138646 | 100121096 | RUSSLEY ROAD | N | RUSSLEY ROAD | 1 | 3071341 | 4 | None | None | sealed | MULTILINESTRING ((1563395.550 5183811.195, 156... |
2633 rows × 11 columns
%%time
sa2["Hdist_motorway"] = sa2.centroid.distance(mway.dissolve().geometry[0])
ax = sa2.plot(column="Hdist_motorway", legend=True)
gpd.clip(mway, sa2.envelope).plot(ax=ax, color="red")
CPU times: user 2.3 s, sys: 130 ms, total: 2.43 s Wall time: 2.28 s
<AxesSubplot:>
rail = gpd.read_file("input/lds-nz-railway-centrelines-topo-150k-SHP.zip")
rail
| t50_fid | name_ascii | macronated | name | status | track_type | rway_use | veh_type | geometry | |
|---|---|---|---|---|---|---|---|---|---|
| 0 | 100120519 | PALMERSTON NORTH - GISBORNE LINE | N | PALMERSTON NORTH - GISBORNE LINE | disused | single | None | None | LINESTRING (2026018.548 5694000.000, 2026005.7... |
| 1 | 7927728 | None | N | None | disused | single | None | None | LINESTRING (1845787.041 5784240.726, 1845749.1... |
| 2 | 2462744 | PALMERSTON NORTH - GISBORNE LINE | N | PALMERSTON NORTH - GISBORNE LINE | disused | single | None | None | LINESTRING (1972000.000 5666908.227, 1971949.1... |
| 3 | 2462838 | TANEATUA BRANCH | Y | TĀNEATUA BRANCH | disused | single | None | None | LINESTRING (1948000.000 5782634.539, 1947964.5... |
| 4 | 6225482 | KINGSTON BRANCH | N | KINGSTON BRANCH | disused | single | None | None | LINESTRING (1261528.423 4961310.584, 1261546.2... |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 527 | 100118890 | STRATFORD - OKAHUKURA LINE | N | STRATFORD - OKAHUKURA LINE | None | single | None | None | LINESTRING (1732000.000 5657175.250, 1732026.9... |
| 528 | 100118892 | NORTH ISLAND MAIN TRUNK | N | NORTH ISLAND MAIN TRUNK | None | single | None | None | LINESTRING (1794254.935 5696679.903, 1794254.8... |
| 529 | 100121048 | STILLWATER - NGAKAWAU LINE | Y | STILLWATER - NGĀKAWAU LINE | None | single | None | None | LINESTRING (1492000.000 5319908.176, 1491999.4... |
| 530 | 100121049 | STILLWATER - NGAKAWAU LINE | Y | STILLWATER - NGĀKAWAU LINE | None | single | None | None | LINESTRING (1484276.055 5375189.782, 1484283.2... |
| 531 | 100121050 | MURUPARA LINE | N | MURUPARA LINE | None | single | None | None | LINESTRING (1932349.458 5766011.454, 1932334.6... |
532 rows × 9 columns
%%time
sa2["Hdist_rail"] = sa2.centroid.distance(rail.dissolve().geometry[0])
ax = sa2.plot(column="Hdist_rail", legend=True)
gpd.clip(rail, sa2.envelope).plot(ax=ax, color="red")
CPU times: user 3.18 s, sys: 170 ms, total: 3.35 s Wall time: 3.21 s
<AxesSubplot:>
skytower = gpd.GeoSeries(skytower, crs="EPSG:2193")
skytower
0 POINT (1757109.809 5920500.841) dtype: geometry
%%time
result = drive(sa2.centroid.to_crs(epsg=4326), skytower.to_crs(epsg=4326))
CPU times: user 314 ms, sys: 28 µs, total: 314 ms Wall time: 2 s
sa2["Mdist_skytower"] = [r[0] for r in result["distances"]]
sa2["Mtime_skytower"] = [r[0] for r in result["durations"]]
sa2.plot(column="Mdist_skytower", legend=True)
<AxesSubplot:>
sa2[sa2.Mdist_skytower < 5000].plot(column="Mdist_skytower", legend=True)
<AxesSubplot:>
%%time
north = drive(sa2.centroid.translate(yoff=-100).to_crs(epsg=4326), skytower.to_crs(epsg=4326))
west = drive(sa2.centroid.translate(xoff=-100).to_crs(epsg=4326), skytower.to_crs(epsg=4326))
east = drive(sa2.centroid.translate(xoff=100).to_crs(epsg=4326), skytower.to_crs(epsg=4326))
south = drive(sa2.centroid.translate(yoff=100).to_crs(epsg=4326), skytower.to_crs(epsg=4326))
CPU times: user 1.07 s, sys: 30 ms, total: 1.1 s Wall time: 7.23 s
sa2["Mdist_skytower_center"] = [r[0] for r in result["distances"]]
sa2["Mdist_skytower_north"] = [r[0] for r in north["distances"]]
sa2["Mdist_skytower_west"] = [r[0] for r in west["distances"]]
sa2["Mdist_skytower_east"] = [r[0] for r in east["distances"]]
sa2["Mdist_skytower_south"] = [r[0] for r in south["distances"]]
sa2[["Mdist_skytower_center",
"Mdist_skytower_north",
"Mdist_skytower_west",
"Mdist_skytower_east",
"Mdist_skytower_south"]].describe()
| Mdist_skytower_center | Mdist_skytower_north | Mdist_skytower_west | Mdist_skytower_east | Mdist_skytower_south | |
|---|---|---|---|---|---|
| count | 556.000000 | 556.000000 | 556.000000 | 556.000000 | 556.000000 |
| mean | 21711.173022 | 21741.161151 | 21706.032554 | 21707.394424 | 21685.279317 |
| std | 15591.755437 | 15569.675209 | 15597.048762 | 15597.416523 | 15609.189601 |
| min | 393.900000 | 487.300000 | 213.100000 | 401.800000 | 286.000000 |
| 25% | 11998.775000 | 12189.525000 | 12101.775000 | 12018.450000 | 12039.750000 |
| 50% | 18263.150000 | 18205.100000 | 18212.050000 | 18094.900000 | 18156.150000 |
| 75% | 26795.150000 | 26853.250000 | 26910.950000 | 26806.200000 | 26840.200000 |
| max | 136328.400000 | 136328.400000 | 136328.400000 | 136328.400000 | 136328.400000 |
sa2["Mdist_skytower"] = sa2[[
"Mdist_skytower_center",
"Mdist_skytower_north",
"Mdist_skytower_west",
"Mdist_skytower_east",
"Mdist_skytower_south"
]].min(axis=1)
(sa2["Mdist_skytower_center"] - sa2["Mdist_skytower"]).describe()
count 556.000000 mean 240.382914 std 433.754953 min 0.000000 25% 85.350000 50% 118.050000 75% 229.500000 max 5626.500000 dtype: float64
sa2[sa2.Mdist_skytower < 5000].plot(column="Mdist_skytower", legend=True)
plt.title("SA2 <5KM driving from the Skytower")
Text(0.5, 1.0, 'SA2 <5KM driving from the Skytower')
sa2["Mtime_skytower_center"] = [r[0] for r in result["durations"]]
sa2["Mtime_skytower_north"] = [r[0] for r in north["durations"]]
sa2["Mtime_skytower_west"] = [r[0] for r in west["durations"]]
sa2["Mtime_skytower_east"] = [r[0] for r in east["durations"]]
sa2["Mtime_skytower_south"] = [r[0] for r in south["durations"]]
display(sa2[[
"Mtime_skytower_center",
"Mtime_skytower_north",
"Mtime_skytower_west",
"Mtime_skytower_east",
"Mtime_skytower_south"
]].describe())
sa2["Mtime_skytower"] = sa2[[
"Mtime_skytower_center",
"Mtime_skytower_north",
"Mtime_skytower_west",
"Mtime_skytower_east",
"Mtime_skytower_south"
]].min(axis=1)
sa2["Mtime_skytower"].describe() # Units are seconds
| Mtime_skytower_center | Mtime_skytower_north | Mtime_skytower_west | Mtime_skytower_east | Mtime_skytower_south | |
|---|---|---|---|---|---|
| count | 556.000000 | 556.000000 | 556.000000 | 556.000000 | 556.000000 |
| mean | 1606.751439 | 1608.667266 | 1605.592266 | 1607.388669 | 1603.354856 |
| std | 3311.208272 | 3310.621902 | 3311.907132 | 3311.035801 | 3311.851192 |
| min | 78.300000 | 92.400000 | 46.600000 | 74.200000 | 50.100000 |
| 25% | 888.275000 | 877.325000 | 880.175000 | 884.750000 | 882.300000 |
| 50% | 1207.200000 | 1200.600000 | 1205.600000 | 1190.000000 | 1195.950000 |
| 75% | 1583.625000 | 1593.575000 | 1591.175000 | 1589.000000 | 1569.025000 |
| max | 70470.600000 | 70470.600000 | 70470.600000 | 70470.600000 | 70470.600000 |
count 556.000000 mean 1584.266547 std 3312.415212 min 46.600000 25% 862.250000 50% 1173.800000 75% 1543.900000 max 70470.600000 Name: Mtime_skytower, dtype: float64
sa2[sa2.Mtime_skytower < 60*20].plot(column="Mtime_skytower", legend=True)
plt.title("SA2s < 20 minutes drive from the Skytower")
Text(0.5, 1.0, 'SA2s < 20 minutes drive from the Skytower')
%%time
result = drive(from_points=sa2.centroid.to_crs(epsg=4326), to_points=transit.geometry.to_crs(epsg=4326))
CPU times: user 273 ms, sys: 22 µs, total: 273 ms Wall time: 1.6 s
min_indices = np.argmin(result["distances"], axis=1)
sa2["Mdist_rapid_transit"] = [result["distances"][i][min_indices[i]] for i in range(len(min_indices))]
sa2["Mtime_rapid_transit"] = [result["durations"][i][min_indices[i]] for i in range(len(min_indices))]
sa2["Mdist_rapid_transit_name"] = list(transit.Location[min_indices])
sa2.Mdist_rapid_transit_name.value_counts()
Constellation Bus Station 41 Hibiscus Coast Bus Station 36 Half Moon Bay 30 Manukau 27 Papakura 23 Pukekohe 23 West Harbour 18 Onehunga 17 Britomart 15 Middlemore 15 Glen Innes 14 Manurewa 14 Mt. Albert 13 Avondale 12 Kingsland 12 Otahuhu 12 Sunnyvale 11 Papatoetoe 11 Henderson 11 Homai 10 Sylvia Park 10 Smales Farm Bus Station 9 Fruitvale Road 9 New Lynn 9 Morningside 9 Glen Eden 9 Ranui 9 Swanson 8 Gulf Harbour 8 Sturges Road 8 Beach Haven 7 Mt Eden 7 Akoranga Bus Station 7 Takanini 7 Grafton 6 Newmarket 6 Ellerslie 6 Matiatia 6 Birkenhead 5 Orakei 5 Meadowbank 4 Panmure 4 Puhinui 4 Pine Harbour 4 Bayswater 4 Parnell 3 Remuera 3 Northcote 3 Devonport 3 Baldwin Avenue 2 Hobsonville 2 Te Papapa 2 Stanley Bay 1 Penrose 1 Te Mahia 1 Name: Mdist_rapid_transit_name, dtype: int64
ax = sa2.to_crs(transit.crs).plot(column="Mdist_rapid_transit", legend=True)
transit.plot(ax=ax, color="red")
<AxesSubplot:>
ax = sa2[sa2.Mdist_rapid_transit < 5000].to_crs(transit.crs).plot(column="Mdist_rapid_transit", legend=True)
transit.plot(ax=ax, color="red")
<AxesSubplot:>
%%time
result = drive(from_points=sa2.centroid.to_crs(epsg=4326), to_points=onramps.centroid.to_crs(epsg=4326))
CPU times: user 369 ms, sys: 197 µs, total: 369 ms Wall time: 2.59 s
min_indices = np.argmin(result["distances"], axis=1)
sa2["Mdist_onramp"] = [result["distances"][i][min_indices[i]] for i in range(len(min_indices))]
sa2["Mtime_onramp"] = [result["durations"][i][min_indices[i]] for i in range(len(min_indices))]
sa2["Mdist_onramp_name"] = list(onramps.DESCRIPTIO.iloc[min_indices])
ax = sa2.plot(column="Mdist_onramp", legend=True)
onramps.plot(ax=ax, color="red")
<AxesSubplot:>
ax = sa2[sa2.Mdist_onramp < 5000].plot(column="Mdist_onramp", legend=True)
onramps.plot(ax=ax, color="red")
<AxesSubplot:>
sa2 = sa2.rename(columns={"LAND_AREA_SQ_KM": "Area"})
%%time
for k,v in tqdm(zones_of_interest.items()):
subset = zones[zones.ZONE_resol.str.startswith(k)] # Just the matched zones (residential, business, rural etc)
print(k, len(subset))
clipped = gpd.clip(sa2, subset) # Clip the local areas to the matched zones
sa2[v] = clipped.area # Store the resulting area. Unit is m²
Residential 22398 Residential - Single House Zone 3766 Residential - Mixed Housing Suburban Zone 9775 Residential - Mixed Housing Urban Zone 5864 Residential - Terrace Housing and Apartment Building Zone 2134 Residential - Large Lot Zone 379 Future Urban Zone 282 Hauraki Gulf Islands 221 Business 3237 Rural 2750 Open Space 6667 CPU times: user 10min 40s, sys: 159 ms, total: 10min 40s Wall time: 10min 39s
cols = list(zones_of_interest.values())
sa2[cols] = sa2[cols].fillna(0)
%%time
sa2["Coast_indicator"] = sa2.geometry.progress_apply(lambda poly: coastline.intersects(poly))
CPU times: user 48.9 s, sys: 49.8 ms, total: 49 s Wall time: 48.9 s
sa2.plot(column="Coast_indicator", legend=True, edgecolor="red")
<AxesSubplot:>
pop = pd.read_csv("input/Individual_part1_totalNZ-wide_format_updated_16-7-20.csv")
pop = pop.set_index("Area_code")
pop
| Area_code_and_description | Area_description | Census_2006_usually_resident_population_count | Census_2013_usually_resident_population_count | Census_2018_usually_resident_population_count | Census_2006_census_night_population_count | Census_2013_census_night_population_count | Census_2018_census_night_population_count | Census_2018_Unit_record_data_source_11_2018_Census_individual_form_CURP | Census_2018_Unit_record_data_source_12_Individuals_on_the_household_listing_only_CURP | ... | Census_2013_Maori_descent_04_Dont_know_CURP | Census_2013_Maori_descent_Total_stated_CURP | Census_2013_Maori_descent_99_Not_elsewhere_included_CURP | Census_2013_Maori_descent_Total_CURP | Census_2018_Maori_descent_01_Māori_descent_CURP | Census_2018_Maori_descent_02_No_Māori_descent_CURP | Census_2018_Maori_descent_04_Dont_know_CURP | Census_2018_Maori_descent_Total_stated_CURP | Census_2018_Maori_descent_99_Not_elsewhere_included_CURP | Census_2018_Maori_descent_Total_CURP | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Area_code | |||||||||||||||||||||
| 7000000 | SA1 7000000 | 7000000 | 165 | 144 | 141 | 177 | 150 | 138 | 114 | 6 | ... | 3 | 132 | 12 | 144 | 135 | 6 | 3 | 141 | 0 | 141 |
| 7000001 | SA1 7000001 | 7000001 | 93 | 105 | 114 | 117 | 108 | 120 | 87 | 3 | ... | 3 | 87 | 15 | 105 | 96 | 18 | 0 | 114 | 0 | 114 |
| 7000002 | SA1 7000002 | 7000002 | 0 | 0 | 0 | 0 | 0 | 0 | C | C | ... | C | C | C | 0 | C | C | C | C | C | 0 |
| 7000003 | SA1 7000003 | 7000003 | 216 | 171 | 225 | 219 | 168 | 222 | 171 | 12 | ... | 3 | 144 | 24 | 171 | 210 | 15 | 0 | 225 | 0 | 225 |
| 7000004 | SA1 7000004 | 7000004 | 90 | 102 | 138 | 90 | 105 | 129 | 117 | 3 | ... | 0 | 75 | 27 | 102 | 102 | 30 | 3 | 138 | 0 | 138 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 16 | 16 Tasman Region | Tasman Region | 44625 | 47157 | 52389 | 48306 | 51258 | 55206 | 46428 | 1443 | ... | 1074 | 43251 | 3903 | 47157 | 5589 | 45471 | 1329 | 52389 | 0 | 52389 |
| 17 | 17 Nelson Region | Nelson Region | 42888 | 46437 | 50880 | 45372 | 48444 | 53082 | 44922 | 1680 | ... | 1173 | 42528 | 3909 | 46437 | 6336 | 43167 | 1377 | 50880 | 0 | 50880 |
| 18 | 18 Marlborough Region | Marlborough Region | 42558 | 43416 | 47340 | 46179 | 46302 | 50562 | 41457 | 1515 | ... | 1104 | 39237 | 4179 | 43416 | 7341 | 38583 | 1416 | 47340 | 0 | 47340 |
| 99 | 99 Area Outside Region | Area Outside Region | 618 | 600 | 669 | 660 | 600 | 822 | 618 | 33 | ... | 24 | 546 | 54 | 600 | 432 | 210 | 27 | 669 | 0 | 669 |
| total | Total NZ (Regional Council) | Total NZ (Regional Council) | 4027947 | 4242048 | 4699755 | 4143282 | 4353198 | 4793358 | 3971892 | 203010 | ... | 87234 | 3821445 | 420603 | 4242048 | 869850 | 3715050 | 114855 | 4699755 | 0 | 4699755 |
32521 rows × 435 columns
sa2["Census2018_population"] = pop.Census_2018_usually_resident_population_count[sa2.Code].tolist()
sa2.plot(column="Census2018_population", legend=True)
<AxesSubplot:>
sa2[["Name", "Census2018_population"]].sort_values("Census2018_population", ascending=False).head(20)
| Name | Census2018_population | |
|---|---|---|
| 1000 | Ormiston South | 5514 |
| 982 | Papatoetoe Central | 5085 |
| 1059 | Pukekohe West | 4980 |
| 706 | Aorere South | 4944 |
| 919 | Point England | 4923 |
| 137 | Weymouth North | 4911 |
| 927 | Mount Wellington Central | 4827 |
| 974 | Papatoetoe North | 4779 |
| 788 | Sunnyvale East | 4773 |
| 128 | Papatoetoe West | 4755 |
| 747 | Birkdale South | 4725 |
| 786 | Glendene North | 4686 |
| 837 | Clover Park North | 4683 |
| 1056 | Pukekohe North West | 4674 |
| 749 | Birkenhead West | 4653 |
| 844 | New Lynn South | 4650 |
| 1011 | Homai West | 4605 |
| 796 | Point Chevalier East | 4602 |
| 869 | Mount Roskill White Swan | 4596 |
| 732 | Massey South | 4554 |
sa2.Census2018_population.plot(kind="hist", bins=100)
<AxesSubplot:ylabel='Frequency'>
dwellings = pd.read_excel("input/2018_census_dwellings_by_SA2.xlsx", sheet_name="Table 1", skiprows=11)
dwellings = dwellings[~pd.isna(dwellings.Total)]
dwellings.index = dwellings["Unnamed: 0"].str.split(" ").str[0]
dwellings = dwellings.replace("..C", np.nan)
dwellings
| Unnamed: 0 | Loss | Zero income | $1-$5,000 | $5,001-$10,000 | $10,001-$15,000 | $15,001-$20,000 | $20,001-$25,000 | $25,001-$30,000 | $30,001-$35,000 | ... | $50,001-$60,000 | $60,001-$70,000 | $70,001-$100,000 | $100,001-$150,000 | $150,001 or More | Total stated | Not stated | Median Household\n Income ($) | Mean Household \nIncome ($) | Total | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Unnamed: 0 | |||||||||||||||||||||
| 100100 | 100100 North Cape | 3.0 | 6.0 | 3.0 | 12.0 | 12.0 | 45.0 | 66.0 | 12.0 | 21.0 | ... | 48.0 | 39.0 | 69.0 | 60.0 | 42.0 | 525.0 | 105.0 | 49400.0 | 65300.0 | 630.0 |
| 100200 | 100200 Rangaunu Harbour | 3.0 | 6.0 | 3.0 | 9.0 | 30.0 | 48.0 | 78.0 | 24.0 | 21.0 | ... | 42.0 | 63.0 | 123.0 | 96.0 | 48.0 | 693.0 | 81.0 | 56400.0 | 67000.0 | 774.0 |
| 100300 | 100300 Inlets Far North District | 0.0 | 0.0 | 0.0 | 3.0 | 0.0 | 3.0 | 6.0 | 0.0 | 0.0 | ... | 0.0 | 0.0 | 3.0 | 0.0 | 3.0 | 30.0 | 30.0 | 24400.0 | 41300.0 | 60.0 |
| 100400 | 100400 Karikari Peninsula | 3.0 | 0.0 | 6.0 | 0.0 | 12.0 | 12.0 | 54.0 | 12.0 | 18.0 | ... | 30.0 | 27.0 | 57.0 | 51.0 | 24.0 | 393.0 | 69.0 | 47200.0 | 62800.0 | 465.0 |
| 100500 | 100500 Tangonge | 3.0 | 3.0 | 0.0 | 0.0 | 6.0 | 12.0 | 30.0 | 9.0 | 6.0 | ... | 24.0 | 24.0 | 54.0 | 69.0 | 36.0 | 330.0 | 57.0 | 67100.0 | 81600.0 | 390.0 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 400013 | 400013 Snares Islands | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.0 |
| 400014 | 400014 Oceanic Antipodes Islands | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.0 |
| 400015 | 400015 Antipodes Islands | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.0 |
| 400016 | 400016 Ross Dependency | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.0 |
| Total | Total New Zealand | 4974.0 | 10152.0 | 15312.0 | 14538.0 | 32031.0 | 59175.0 | 105777.0 | 45423.0 | 45849.0 | ... | 107058.0 | 97539.0 | 242262.0 | 295119.0 | 272208.0 | 1526961.0 | 126831.0 | 75700.0 | 91300.0 | 1653792.0 |
2254 rows × 22 columns
sa2["Census2018_dwellings"] = dwellings.Total[sa2.Code].tolist()
sa2["Census2018_avg_HH_income"] = dwellings["Mean Household \nIncome ($)"][sa2.Code].tolist()
sa2["Census2018_med_HH_income"] = dwellings["Median Household\n Income ($)"][sa2.Code].tolist()
sa2.plot(column="Census2018_dwellings", legend=True)
sa2.plot(column="Census2018_avg_HH_income", legend=True)
sa2.plot(column="Census2018_med_HH_income", legend=True)
<AxesSubplot:>
sa2.Census2018_med_HH_income.plot(kind="hist", bins=100)
<AxesSubplot:ylabel='Frequency'>
sa2
| Code | Name | REGC2018_V1_00 | REGC2018_V1_00_NAME | TA2018_V1_00 | TA2018_V1_00_NAME | Area | AREA_SQ_KM | Shape_Length | geometry | ... | FU_area | HGI_area | Business_area | Rural_area | Open_area | Coast_indicator | Census2018_population | Census2018_dwellings | Census2018_avg_HH_income | Census2018_med_HH_income | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 33 | 110900 | Dome Valley-Matakana | 02 | Auckland Region | 076 | Auckland | 84.749966 | 84.749966 | 70031.926112 | MULTIPOLYGON (((1749212.212 5981874.734, 17492... | ... | 0.000000e+00 | 0.0 | 7.881657e+04 | 7.396202e+07 | 7.310753e+06 | True | 1530 | 546.0 | 104300.0 | 88200.0 |
| 34 | 111100 | Warkworth West | 02 | Auckland Region | 076 | Auckland | 11.359752 | 11.359752 | 22629.664017 | MULTIPOLYGON (((1748589.106 5974509.939, 17485... | ... | 3.210818e+06 | 0.0 | 1.123969e+06 | 4.675043e+06 | 5.981430e+05 | False | 2481 | 960.0 | 82000.0 | 59600.0 |
| 35 | 111200 | Puhoi Valley | 02 | Auckland Region | 076 | Auckland | 236.529965 | 236.529965 | 142134.544465 | MULTIPOLYGON (((1740248.959 5974374.891, 17402... | ... | 4.191487e+06 | 0.0 | 5.249019e+03 | 2.082993e+08 | 1.175816e+07 | True | 3702 | 1242.0 | 112300.0 | 104100.0 |
| 36 | 112600 | Waikoukou Valley | 02 | Auckland Region | 076 | Auckland | 52.854583 | 52.854583 | 42565.729322 | MULTIPOLYGON (((1735693.112 5942203.564, 17359... | ... | 0.000000e+00 | 0.0 | 0.000000e+00 | 5.134617e+07 | 4.289403e+05 | False | 1728 | 555.0 | 132100.0 | 126700.0 |
| 37 | 112800 | Hatfields Beach | 02 | Auckland Region | 076 | Auckland | 0.923255 | 0.923255 | 4314.440531 | MULTIPOLYGON (((1750981.811 5952164.661, 17512... | ... | 5.546933e+04 | 0.0 | 0.000000e+00 | 2.666470e+00 | 2.052774e+05 | True | 1554 | 558.0 | 108600.0 | 97100.0 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 1060 | 165600 | Cape Hill | 02 | Auckland Region | 076 | Auckland | 0.717375 | 0.717375 | 3569.522355 | MULTIPOLYGON (((1770139.188 5882466.221, 17701... | ... | 1.126037e+05 | 0.0 | 0.000000e+00 | 0.000000e+00 | 8.655591e+04 | False | 1551 | 480.0 | 123400.0 | 116300.0 |
| 1061 | 165800 | Rooseville Park | 02 | Auckland Region | 076 | Auckland | 1.345809 | 1.345809 | 6010.945714 | MULTIPOLYGON (((1770089.518 5882192.248, 17701... | ... | 0.000000e+00 | 0.0 | 3.491631e+04 | 0.000000e+00 | 2.654751e+05 | False | 2562 | 864.0 | 104600.0 | 96600.0 |
| 1062 | 166000 | Pukekohe Central | 02 | Auckland Region | 076 | Auckland | 2.939129 | 2.939129 | 9211.796424 | MULTIPOLYGON (((1769469.158 5881504.407, 17694... | ... | 2.698912e+05 | 0.0 | 9.618075e+05 | 0.000000e+00 | 2.140046e+04 | False | 369 | 153.0 | 64500.0 | 53200.0 |
| 1063 | 166300 | Bombay Hills | 02 | Auckland Region | 076 | Auckland | 31.682020 | 31.682020 | 37482.345738 | MULTIPOLYGON (((1777896.083 5888522.045, 17779... | ... | 0.000000e+00 | 0.0 | 8.384705e+04 | 2.835952e+07 | 3.042428e+05 | False | 1974 | 669.0 | 122500.0 | 116200.0 |
| 1064 | 166400 | Ararimu | 02 | Auckland Region | 076 | Auckland | 91.308241 | 91.308241 | 57552.107386 | MULTIPOLYGON (((1786544.313 5891035.131, 17865... | ... | 0.000000e+00 | 0.0 | 1.575686e+03 | 8.638718e+07 | 2.033428e+03 | False | 2124 | 675.0 | 130100.0 | 122600.0 |
556 rows × 50 columns
sa2.drop(columns="geometry").to_csv("output/2018_Statistical_Area_2.csv")
au = gpd.read_file("input/statsnzarea-unit-2013-FGDB.zip!area-unit-2013.gdb")
au = au.rename(columns={"AU2013_V1_00_NAME": "Name", "AU2013_V1_00": "Code"})
au["Centroid_lon"] = au.centroid.to_crs(epsg=4326).x
au["Centroid_lat"] = au.centroid.to_crs(epsg=4326).y
au["Centroid_local_area"] = au.representative_point().apply(lambda point: df.Name[df.contains(point)].iloc[0] if any(df.contains(point)) else np.nan)
problems = au.Code.isin(["520602", "521135", "506634"])
au.at[problems, "Centroid_local_area"] = au[problems].geometry.apply(lambda poly: df.Name[df.intersects(poly)].iloc[0])
au["Centroid_local_area"].value_counts()
Henderson-Massey 29 Manurewa 21 Orakei 20 Auckland Central 20 Mangere-Otahuhu 19 Beach Haven-Birkenhead-Northcote 19 Hibiscus Coast 17 Whau 16 Rodney-Helensville 15 Upper Harbour Local Board Area 15 Papakura 14 Devonport-Takapuna 13 Mt Eden 13 Botany 13 Franklin-Pukekohe 12 Titirangi 12 Rodney-Warkworth 12 Papatoetoe 11 Pakuranga 11 Puketapapa 11 Tamaki 11 Mt Albert 11 East Coast Bays 11 Franklin-Beachlands-Hunua 10 Howick 10 Maungakiekie 9 Otara 9 Waitakere 8 Franklin-Waiuku 5 Rodney-Kumeu-Riverhead 4 Waiheke 3 Rodney-Wellsford 3 Rodney-Dairy Flat 2 Name: Centroid_local_area, dtype: int64
print(f'{sum(~pd.isna(au["Centroid_local_area"]))} matches, {sum(pd.isna(au["Centroid_local_area"]))} unmatched')
409 matches, 1595 unmatched
au = au[~pd.isna(au.Centroid_local_area)]
ax = au.plot(column="Centroid_local_area", cmap="Spectral", legend=True)
au["Hdist_skytower"] = au.distance(skytower.iloc[0])
au.plot(column="Hdist_skytower", legend=True)
<AxesSubplot:>
%%time
au["Hdist_motorway"] = au.centroid.distance(mway.dissolve().geometry[0])
ax = au.plot(column="Hdist_motorway", legend=True)
gpd.clip(mway, au.envelope).plot(ax=ax, color="red")
CPU times: user 1.84 s, sys: 281 ms, total: 2.13 s Wall time: 1.98 s
<AxesSubplot:>
%%time
au["Hdist_rail"] = au.centroid.distance(rail.dissolve().geometry[0])
ax = au.plot(column="Hdist_rail", legend=True)
gpd.clip(rail, au.envelope).plot(ax=ax, color="red")
CPU times: user 2.92 s, sys: 150 ms, total: 3.07 s Wall time: 2.92 s
<AxesSubplot:>
%%time
result = drive(au.centroid.to_crs(epsg=4326), skytower.to_crs(epsg=4326))
north = drive(au.centroid.translate(yoff=-100).to_crs(epsg=4326), skytower.to_crs(epsg=4326))
west = drive(au.centroid.translate(xoff=-100).to_crs(epsg=4326), skytower.to_crs(epsg=4326))
east = drive(au.centroid.translate(xoff=100).to_crs(epsg=4326), skytower.to_crs(epsg=4326))
south = drive(au.centroid.translate(yoff=100).to_crs(epsg=4326), skytower.to_crs(epsg=4326))
CPU times: user 1.16 s, sys: 9.82 ms, total: 1.17 s Wall time: 5.85 s
au["Mdist_skytower_center"] = [r[0] for r in result["distances"]]
au["Mdist_skytower_north"] = [r[0] for r in north["distances"]]
au["Mdist_skytower_west"] = [r[0] for r in west["distances"]]
au["Mdist_skytower_east"] = [r[0] for r in east["distances"]]
au["Mdist_skytower_south"] = [r[0] for r in south["distances"]]
display(au[[
"Mdist_skytower_center",
"Mdist_skytower_north",
"Mdist_skytower_west",
"Mdist_skytower_east",
"Mdist_skytower_south"
]].describe())
au["Mdist_skytower"] = au[[
"Mdist_skytower_center",
"Mdist_skytower_north",
"Mdist_skytower_west",
"Mdist_skytower_east",
"Mdist_skytower_south"
]].min(axis=1)
display(au["Mdist_skytower"].describe())
au["Mtime_skytower_center"] = [r[0] for r in result["durations"]]
au["Mtime_skytower_north"] = [r[0] for r in north["durations"]]
au["Mtime_skytower_west"] = [r[0] for r in west["durations"]]
au["Mtime_skytower_east"] = [r[0] for r in east["durations"]]
au["Mtime_skytower_south"] = [r[0] for r in south["durations"]]
display(au[[
"Mtime_skytower_center",
"Mtime_skytower_north",
"Mtime_skytower_west",
"Mtime_skytower_east",
"Mtime_skytower_south"
]].describe())
au["Mtime_skytower"] = au[[
"Mtime_skytower_center",
"Mtime_skytower_north",
"Mtime_skytower_west",
"Mtime_skytower_east",
"Mtime_skytower_south"
]].min(axis=1)
display(au["Mtime_skytower"].describe())
| Mdist_skytower_center | Mdist_skytower_north | Mdist_skytower_west | Mdist_skytower_east | Mdist_skytower_south | |
|---|---|---|---|---|---|
| count | 409.000000 | 409.000000 | 409.000000 | 409.000000 | 409.000000 |
| mean | 22795.600244 | 22781.618093 | 22788.159902 | 22784.258191 | 22748.530073 |
| std | 16379.256603 | 16359.032083 | 16357.011292 | 16376.317853 | 16403.724191 |
| min | 557.400000 | 715.700000 | 707.700000 | 1143.700000 | 446.000000 |
| 25% | 12472.900000 | 12402.700000 | 12501.100000 | 12532.200000 | 12441.400000 |
| 50% | 18613.100000 | 18513.300000 | 18539.800000 | 18517.400000 | 18512.800000 |
| 75% | 27316.300000 | 27355.900000 | 27221.000000 | 27385.600000 | 27368.200000 |
| max | 89091.100000 | 89018.800000 | 89160.200000 | 87888.800000 | 89163.300000 |
count 409.000000 mean 22531.981907 std 16385.776202 min 446.000000 25% 12227.900000 50% 18255.700000 75% 27187.300000 max 87888.800000 Name: Mdist_skytower, dtype: float64
| Mtime_skytower_center | Mtime_skytower_north | Mtime_skytower_west | Mtime_skytower_east | Mtime_skytower_south | |
|---|---|---|---|---|---|
| count | 409.000000 | 409.000000 | 409.000000 | 409.000000 | 409.000000 |
| mean | 1483.669438 | 1484.317115 | 1483.113692 | 1483.419804 | 1481.469193 |
| std | 1264.627438 | 1263.705290 | 1262.503954 | 1264.549168 | 1264.771003 |
| min | 99.900000 | 110.600000 | 148.500000 | 172.300000 | 73.100000 |
| 25% | 910.500000 | 926.500000 | 910.300000 | 913.200000 | 905.600000 |
| 50% | 1225.300000 | 1222.200000 | 1228.600000 | 1219.700000 | 1226.300000 |
| 75% | 1621.600000 | 1623.600000 | 1603.700000 | 1625.500000 | 1615.500000 |
| max | 14754.400000 | 14754.400000 | 14753.900000 | 14802.600000 | 14754.400000 |
count 409.000000 mean 1459.823227 std 1263.493685 min 73.100000 25% 883.600000 50% 1203.900000 75% 1588.300000 max 14753.900000 Name: Mtime_skytower, dtype: float64
au.plot(column="Mdist_skytower", legend=True)
<AxesSubplot:>
au[au.Mdist_skytower < 10000].plot(column="Mdist_skytower", legend=True)
<AxesSubplot:>
%%time
result = drive(from_points=au.centroid.to_crs(epsg=4326), to_points=transit.geometry.to_crs(epsg=4326))
CPU times: user 212 ms, sys: 39.6 ms, total: 252 ms Wall time: 1.24 s
min_indices = np.argmin(result["distances"], axis=1)
au["Mdist_rapid_transit"] = [result["distances"][i][min_indices[i]] for i in range(len(min_indices))]
au["Mtime_rapid_transit"] = [result["durations"][i][min_indices[i]] for i in range(len(min_indices))]
au["Mdist_rapid_transit_name"] = list(transit.Location[min_indices])
au.Mdist_rapid_transit_name.value_counts()
Hibiscus Coast Bus Station 36 Constellation Bus Station 25 Half Moon Bay 24 Papakura 21 West Harbour 18 Pukekohe 15 Manukau 13 Papatoetoe 12 Otahuhu 12 Kingsland 11 Manurewa 11 Ranui 10 Onehunga 10 Glen Innes 10 Sunnyvale 10 Mt. Albert 10 Henderson 9 Middlemore 9 Glen Eden 9 Swanson 8 Avondale 7 Morningside 7 Sylvia Park 6 Homai 6 New Lynn 6 Beach Haven 6 Fruitvale Road 6 Grafton 6 Britomart 5 Smales Farm Bus Station 5 Akoranga Bus Station 5 Orakei 5 Newmarket 4 Gulf Harbour 4 Birkenhead 4 Ellerslie 4 Mt Eden 4 Panmure 3 Te Papapa 3 Remuera 3 Bayswater 3 Sturges Road 3 Takanini 3 Parnell 2 Matiatia 2 Devonport 2 Northcote 2 Puhinui 2 Pine Harbour 2 Baldwin Avenue 2 Stanley Bay 1 Penrose 1 Meadowbank 1 Te Mahia 1 Name: Mdist_rapid_transit_name, dtype: int64
ax = au.to_crs(transit.crs).plot(column="Mdist_rapid_transit", legend=True)
transit.plot(ax=ax, color="red")
<AxesSubplot:>
ax = au[au.Mdist_rapid_transit < 10000].to_crs(transit.crs).plot(column="Mdist_rapid_transit", legend=True)
transit.plot(ax=ax, color="red")
<AxesSubplot:>
%%time
result = drive(from_points=au.centroid.to_crs(epsg=4326), to_points=onramps.centroid.to_crs(epsg=4326))
CPU times: user 313 ms, sys: 29.6 ms, total: 343 ms Wall time: 1.82 s
min_indices = np.argmin(result["distances"], axis=1)
au["Mdist_onramp"] = [result["distances"][i][min_indices[i]] for i in range(len(min_indices))]
au["Mtime_onramp"] = [result["durations"][i][min_indices[i]] for i in range(len(min_indices))]
au["Mdist_onramp_name"] = list(onramps.DESCRIPTIO.iloc[min_indices])
ax = au.plot(column="Mdist_onramp", legend=True)
onramps.plot(ax=ax, color="red")
<AxesSubplot:>
ax = au[au.Mdist_onramp < 5000].plot(column="Mdist_onramp", legend=True)
onramps.plot(ax=ax, color="red")
<AxesSubplot:>
au = au.rename(columns={"LAND_AREA_SQ_KM": "Area"})
%%time
for k,v in tqdm(zones_of_interest.items()):
subset = zones[zones.ZONE_resol.str.startswith(k)] # Just the matched zones (residential, business, rural etc)
print(k, len(subset))
clipped = gpd.clip(au, subset) # Clip the local areas to the matched zones
au[v] = clipped.area # Store the resulting area. Unit is m²
Residential 22398 Residential - Single House Zone 3766 Residential - Mixed Housing Suburban Zone 9775 Residential - Mixed Housing Urban Zone 5864 Residential - Terrace Housing and Apartment Building Zone 2134 Residential - Large Lot Zone 379 Future Urban Zone 282 Hauraki Gulf Islands 221 Business 3237 Rural 2750 Open Space 6667 CPU times: user 10min 2s, sys: 1.97 s, total: 10min 3s Wall time: 10min 3s
au = au.fillna(0)
%%time
au["Coast_indicator"] = au.geometry.progress_apply(lambda poly: coastline.intersects(poly))
CPU times: user 30.7 s, sys: 150 ms, total: 30.8 s Wall time: 30.7 s
au.plot(column="Coast_indicator", legend=True, edgecolor="red")
<AxesSubplot:>
pop = pd.read_csv("input/2013-mb-dataset-Total-New-Zealand-Individual-Part-1.csv", encoding='unicode_escape', low_memory=False)
pop
| Area_Code_and_Description | Code | Description | 2001_Census_census_usually_resident_population_count(1) | 2006_Census_census_usually_resident_population_count(1) | 2013_Census_census_usually_resident_population_count(1) | 2001_Census_census_night_population_count(2) | 2006_Census_census_night_population_count(2) | 2013_Census_census_night_population_count(2) | 2001_Census_sex_for_the_census_usually_resident_population_count(1)_Male | ... | 2006_Census_Maori_descent_for_the_census_usually_resident_population_count(1)_Don't_Know | 2006_Census_Maori_descent_for_the_census_usually_resident_population_count(1)_Total_people_stated | 2006_Census_Maori_descent_for_the_census_usually_resident_population_count(1)_Not_Elsewhere_Included(14) | 2006_Census_Maori_descent_for_the_census_usually_resident_population_count(1)_Total_people | 2013_Census_Maori_descent_for_the_census_usually_resident_population_count(1)_Maori_Descent | 2013_Census_Maori_descent_for_the_census_usually_resident_population_count(1)_No_Maori_Descent | 2013_Census_Maori_descent_for_the_census_usually_resident_population_count(1)_Don't_Know | 2013_Census_Maori_descent_for_the_census_usually_resident_population_count(1)_Total_people_stated | 2013_Census_Maori_descent_for_the_census_usually_resident_population_count(1)_Not_Elsewhere_Included(14) | 2013_Census_Maori_descent_for_the_census_usually_resident_population_count(1)_Total_people | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | MB 0000100 | 0000100 | NaN | 9.0 | 3.0 | 3.0 | 15.0 | 9.0 | 9.0 | 6 | ... | ..C | ..C | ..C | 3.0 | ..C | ..C | ..C | ..C | ..C | 3.0 |
| 1 | MB 0000200 | 0000200 | NaN | 90.0 | 87.0 | 84.0 | 99.0 | 96.0 | 84.0 | 42 | ... | 0 | 81 | 6 | 90.0 | 69 | 3 | 3 | 75 | 9 | 84.0 |
| 2 | MB 0000300 | 0000300 | NaN | 90.0 | 72.0 | 57.0 | 93.0 | 72.0 | 54.0 | 42 | ... | 0 | 66 | 6 | 75.0 | 54 | 0 | 0 | 57 | 3 | 57.0 |
| 3 | MB 0000400 | 0000400 | NaN | 30.0 | 21.0 | 30.0 | 45.0 | 39.0 | 30.0 | 18 | ... | 0 | 21 | 3 | 21.0 | 12 | 6 | 3 | 21 | 6 | 27.0 |
| 4 | MB 0000501 | 0000501 | NaN | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | ..C | ... | ..C | ..C | ..C | 0.0 | ..C | ..C | ..C | ..C | ..C | 0.0 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 49028 | Footnotes for specific variables or categories | 13 | Consists of people who were too young to talk ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
| 49029 | Footnotes for specific variables or categories | 14 | Consists of response unidentifiable and not st... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
| 49030 | Symbols | ..C | confidential | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
| 49031 | Symbols | * | not able to be calculated | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
| 49032 | Source | Statistics New Zealand | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
49033 rows × 227 columns
pop = pop.set_index("Code")
au["Census2013_population"] = pop["2013_Census_census_usually_resident_population_count(1)"][au.Code].tolist()
au.plot(column="Census2013_population", legend=True)
<AxesSubplot:>
households = households.set_index("Code")
au["Census2013_dwellings"] = households["2013_Census_total_households_in_occupied_private_dwellings"][au.Code].tolist()
au.plot(column="Census2013_dwellings", legend=True)
<AxesSubplot:>
# Can't do this one due to missing data
key
'2013_Census_total_household_income_(grouped)(2)(3)(4)_for_households_in_occupied_private_dwellings_Median_household_income_($)(18)(23)'
au["Census2013_med_HH_income"] = households[key][au.Code].replace("..C", np.nan).astype(float).tolist()
au.plot(column="Census2013_med_HH_income", legend=True)
<AxesSubplot:>
au
| Code | Name | AREA_SQ_KM | Area | Shape_Length | geometry | Centroid_lon | Centroid_lat | Centroid_local_area | Hdist_skytower | ... | LL_area | FU_area | HGI_area | Business_area | Rural_area | Open_area | Coast_indicator | Census2013_population | Census2013_dwellings | Census2013_med_HH_income | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 30 | 505604 | Kumeu East | 13.475480 | 13.310426 | 24949.702400 | MULTIPOLYGON (((1741052.061 5931333.839, 17409... | 174.575870 | -36.770909 | Rodney-Kumeu-Riverhead | 16136.289703 | ... | 0.000000e+00 | 1.548950e+06 | 0.00000 | 142889.422393 | 1.037018e+07 | 216410.264282 | True | 1449.0 | 474.0 | 87700.0 |
| 31 | 505605 | Kumeu West | 24.388271 | 24.388271 | 32928.601801 | MULTIPOLYGON (((1734788.301 5932715.055, 17348... | 174.524933 | -36.774538 | Rodney-Kumeu-Riverhead | 18120.060194 | ... | 0.000000e+00 | 5.636526e+06 | 0.00000 | 0.000000 | 1.749381e+07 | 228438.588353 | False | 1725.0 | 609.0 | 95300.0 |
| 32 | 505804 | Hatfields Beach | 0.908867 | 0.908867 | 5395.241593 | MULTIPOLYGON (((1751816.238 5952425.382, 17518... | 174.690104 | -36.569595 | Hibiscus Coast | 30962.434695 | ... | 5.866738e+04 | 5.548376e+04 | 0.00000 | 0.000000 | 3.031425e+00 | 239098.186254 | True | 1380.0 | 492.0 | 74100.0 |
| 33 | 505805 | Orewa | 5.572505 | 4.229676 | 19755.382361 | MULTIPOLYGON (((1750518.569 5950940.087, 17505... | 174.686690 | -36.589013 | Hibiscus Coast | 27831.510696 | ... | 0.000000e+00 | 5.580409e+01 | 0.00000 | 198410.958707 | 1.820051e+01 | 667506.890608 | True | 8520.0 | 3843.0 | 47600.0 |
| 34 | 505902 | Manly | 3.188420 | 3.188420 | 12908.971485 | MULTIPOLYGON (((1756245.400 5945964.220, 17562... | 174.759445 | -36.630419 | Hibiscus Coast | 22916.935706 | ... | 1.725715e+04 | 0.000000e+00 | 0.00000 | 9822.350907 | 0.000000e+00 | 373513.517145 | True | 6552.0 | 2535.0 | 63800.0 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 967 | 525922 | Bledisloe Park | 9.643551 | 9.643551 | 16858.446790 | MULTIPOLYGON (((1770838.870 5881130.225, 17709... | 174.915328 | -37.215688 | Franklin-Pukekohe | 41041.238084 | ... | 0.000000e+00 | 3.229772e+06 | 0.00000 | 796129.574379 | 1.976364e+06 | 280880.995489 | False | 2934.0 | 1092.0 | 70100.0 |
| 968 | 526103 | Waiuku West | 2.198567 | 2.198567 | 9682.515485 | MULTIPOLYGON (((1753228.371 5877597.817, 17532... | 174.724591 | -37.248157 | Franklin-Waiuku | 42771.770108 | ... | 2.023937e+04 | 0.000000e+00 | 0.00000 | 208423.241402 | 2.600376e+05 | 265772.078309 | True | 3102.0 | 1221.0 | 56900.0 |
| 969 | 526104 | Waiuku East | 3.200239 | 3.158494 | 9009.196620 | MULTIPOLYGON (((1753749.548 5877735.141, 17537... | 174.738476 | -37.247818 | Franklin-Waiuku | 42896.783550 | ... | 7.459611e+05 | 0.000000e+00 | 0.00000 | 175045.781248 | 0.000000e+00 | 531138.555183 | True | 3498.0 | 1248.0 | 60500.0 |
| 970 | 526105 | South Waiuku | 2.247383 | 2.247383 | 9084.157884 | MULTIPOLYGON (((1755346.188 5875198.493, 17551... | 174.739741 | -37.259971 | Franklin-Waiuku | 44558.447057 | ... | 1.185644e+06 | 0.000000e+00 | 0.00000 | 0.000000 | 5.934895e+05 | 7326.469987 | False | 1599.0 | 546.0 | 85600.0 |
| 1962 | 616300 | Browns Island | 0.595336 | 0.595336 | 3430.690993 | MULTIPOLYGON (((1769155.864 5922005.512, 17691... | 174.894441 | -36.830617 | Waiheke | 11500.202550 | ... | 0.000000e+00 | 0.000000e+00 | 586560.33176 | 0.000000 | 0.000000e+00 | 0.000000 | True | 0.0 | 0.0 | NaN |
409 rows × 45 columns
au.drop(columns="geometry").to_csv("output/2013_Area_Unit.csv")
zones
| OBJECTID | CONTOUR | created_da | DocumentUR | GlobalID | GROUPZONE | GROUPZONE_ | last_edite | NAME | PARCEL_BAS | ... | TYPE | TYPE_resol | VERSIONSTA | VERSIONS_1 | ZONE | ZONE_resol | ZONEHEIGHT | SHAPE_Leng | SHAPE_Area | geometry | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 1.0 | None | 20160718211 | None | {4C8F9436-7EA6-417E-B64F-15FCD44459F6} | 2 | Residential | 20161111010 | None | None | ... | None | None | 4 | Operative | 60 | Residential - Mixed Housing Urban Zone | NaN | 285.664016 | 2.050275e+03 | POLYGON ((1768030.306 5901206.846, 1768033.070... |
| 1 | 2.0 | None | 20160718211 | None | {604AAD87-8ED4-4111-8276-47CEE7E81F92} | 1 | Public Open Space | 20161111010 | None | None | ... | None | None | 4 | Operative | 33 | Open Space - Sport and Active Recreation Zone | NaN | 1246.837757 | 1.684599e+04 | POLYGON ((1764267.286 5919989.370, 1764218.153... |
| 2 | 3.0 | None | 20160718211 | None | {8D827DA8-BC5B-437A-B17A-532354F7D037} | 4 | Rural | 20161111010 | None | None | ... | None | None | 4 | Operative | 11 | Rural - Mixed Rural Zone | NaN | 3582.113246 | 6.841744e+05 | POLYGON ((1740091.195 5928308.839, 1740089.844... |
| 3 | 4.0 | None | 20160718211 | None | {96C9E266-3341-4C71-94F1-325F2EE45732} | 2 | Residential | 20161111010 | None | None | ... | None | None | 4 | Operative | 23 | Residential - Large Lot Zone | NaN | 317.098469 | 6.024226e+03 | POLYGON ((1750502.125 5928677.635, 1750456.131... |
| 4 | 5.0 | None | 20160718211 | None | {90B50FEE-45A3-4E88-819A-370751ACDE3D} | 1 | Public Open Space | 20161111010 | None | None | ... | None | None | 4 | Operative | 31 | Open Space - Conservation Zone | NaN | 230.836636 | 5.639794e+02 | POLYGON ((1741460.217 5918191.600, 1741476.745... |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 130295 | 130296.0 | None | 20161115151 | None | {B2F0FB45-80F6-41ED-AA57-A914B195B31E} | 1 | Public Open Space | 20161115151 | None | None | ... | None | None | 4 | Operative | 31 | Open Space - Conservation Zone | NaN | 179382.307787 | 9.040193e+07 | POLYGON ((1738072.631 5902593.421, 1738095.386... |
| 130296 | 130297.0 | None | 20161115151 | None | {A3F4EF61-9162-43C3-90BD-DC84D93C6A64} | 2 | Residential | 20161115151 | None | None | ... | None | None | 4 | Operative | 18 | Residential - Mixed Housing Suburban Zone | NaN | 2769.292040 | 2.920602e+05 | POLYGON ((1770189.642 5905789.775, 1770199.407... |
| 130297 | 130298.0 | None | 20161115151 | None | {A5FC7EB5-76E5-414C-96DD-715C671764B4} | 4 | Rural | 20161115151 | Kaipara South Head and Harbour coastal area | None | ... | None | None | 4 | Operative | 46 | Rural - Rural Coastal Zone | NaN | 27168.762393 | 6.344938e+06 | POLYGON ((1730253.401 5955767.484, 1730254.399... |
| 130298 | 130299.0 | None | 20161115151 | None | {C371B83C-35CE-46A1-B94F-F1E592F271F7} | 2 | Residential | 20161115151 | None | None | ... | None | None | 4 | Operative | 8 | Residential - Terrace Housing and Apartment Bu... | NaN | 1812.109533 | 1.536829e+05 | POLYGON ((1743032.759 5924130.564, 1742996.886... |
| 130299 | 130300.0 | None | 20161115151 | None | {31281924-C74D-4038-8260-FE6D886F4C94} | 2 | Residential | 20161115151 | None | None | ... | None | None | 4 | Operative | 18 | Residential - Mixed Housing Suburban Zone | NaN | 2511.267945 | 3.539606e+05 | POLYGON ((1748979.767 5925771.895, 1748986.615... |
130300 rows × 31 columns
cols = [
"Residential - Terrace Housing and Apartment Building Zone",
"Residential - Mixed Housing Urban Zone",
"Residential - Mixed Housing Suburban Zone"
]
coastline = gpd.read_file("input/lds-nz-coastlines-and-islands-polygons-topo-150k-FGDB.zip!nz-coastlines-and-islands-polygons-topo-150k.gdb")
cmap = matplotlib.colors.ListedColormap(['red', 'blue', 'green'])
ax = zones[zones.ZONE_resol.isin(cols)].plot(column="ZONE_resol", legend=True, cmap=cmap)
xmin, xmax = ax.get_xlim()
ymin, ymax = ax.get_ylim()
print(xmin, xmax, ymin, ymax)
coastline.boundary.plot(ax=ax, color="gray", linewidth=.5)
ax.set_xlim(xmin, xmax)
ax.set_ylim(ymin + 1e4, ymax - 4e4)
ax.add_artist(ScaleBar(1, location="lower left"))
x, y, arrow_length = 0.05, .98, 0.05
ax.annotate('N', xy=(x, y), xytext=(x, y-arrow_length),
arrowprops=dict(facecolor='black', width=5, headwidth=15),
ha='center', va='center', fontsize=20,
xycoords=ax.transAxes)
ax.set_xticks([])
ax.set_yticks([])
1735024.9733500003 1779731.6056499996 5869626.39383 5980785.90757
[]
ax = zones[zones.ZONE_resol == "Residential - Single House Zone"].plot(column="ZONE_resol", legend=True, cmap=cmap)
xmin, xmax = ax.get_xlim()
ymin, ymax = ax.get_ylim()
print(xmin, xmax, ymin, ymax)
coastline.boundary.plot(ax=ax, color="gray", linewidth=.5)
ax.set_xlim(xmin + 1e4, xmax - 1e4)
ax.set_ylim(ymin + 1e4, ymax - 4e4)
ax.add_artist(ScaleBar(1, location="lower left"))
x, y, arrow_length = 0.05, .98, 0.05
ax.annotate('N', xy=(x, y), xytext=(x, y-arrow_length),
arrowprops=dict(facecolor='black', width=5, headwidth=15),
ha='center', va='center', fontsize=20,
xycoords=ax.transAxes)
ax.set_xticks([])
ax.set_yticks([])
1724205.46675 1797169.8142499996 5869706.9597 5992439.5543
[]
df.Name[df.Hdist_skytower < 2e4].tolist()
['Auckland Central', 'Beach Haven-Birkenhead-Northcote', 'Botany', 'Devonport-Takapuna', 'East Coast Bays', 'Henderson-Massey', 'Howick', 'Mangere-Otahuhu', 'Maungakiekie', 'Mt Albert', 'Mt Eden', 'Orakei', 'Otara', 'Pakuranga', 'Papatoetoe', 'Puketapapa', 'Tamaki', 'Titirangi', 'Upper Harbour Local Board Area', 'Whau']
df["Core"] = df.Name.isin(['Auckland Central',
'Beach Haven-Birkenhead-Northcote',
'Botany',
'Devonport-Takapuna',
'East Coast Bays',
'Henderson-Massey',
'Howick',
'Mangere-Otahuhu',
'Manurewa',
'Maungakiekie',
'Mt Albert',
'Mt Eden',
'Orakei',
'Otara',
'Pakuranga',
'Papatoetoe',
'Papakura',
'Puketapapa',
'Tamaki',
'Titirangi',
'Upper Harbour Local Board Area',
'Whau'])
display(df[["Name", "Core"]])
| Name | Core | |
|---|---|---|
| OBJECTID | ||
| 1 | Auckland Central | True |
| 2 | Beach Haven-Birkenhead-Northcote | True |
| 3 | Botany | True |
| 4 | Devonport-Takapuna | True |
| 5 | East Coast Bays | True |
| 6 | Franklin-Beachlands-Hunua | False |
| 7 | Franklin-Pukekohe | False |
| 8 | Franklin-Waiuku | False |
| 9 | Henderson-Massey | True |
| 10 | Hibiscus Coast | False |
| 11 | Howick | True |
| 12 | Mangere-Otahuhu | True |
| 13 | Manurewa | True |
| 14 | Maungakiekie | True |
| 15 | Mt Albert | True |
| 16 | Mt Eden | True |
| 17 | Orakei | True |
| 18 | Otara | True |
| 19 | Pakuranga | True |
| 20 | Papakura | True |
| 21 | Papatoetoe | True |
| 22 | Puketapapa | True |
| 23 | Rodney-Dairy Flat | False |
| 24 | Rodney-Helensville | False |
| 25 | Rodney-Kumeu-Riverhead | False |
| 26 | Rodney-Warkworth | False |
| 27 | Rodney-Wellsford | False |
| 28 | Tamaki | True |
| 29 | Titirangi | True |
| 30 | Upper Harbour Local Board Area | True |
| 31 | Waiheke | False |
| 32 | Waitakere | False |
| 33 | Whau | True |
cmap = matplotlib.colors.ListedColormap(['orange', 'yellow'])
ax = df.plot(column = "Core", cmap=cmap, edgecolor="black", linewidth=.5)
xmin, xmax = ax.get_xlim()
ymin, ymax = ax.get_ylim()
print(xmin, xmax, ymin, ymax)
coastline.boundary.plot(ax=ax, color="gray", linewidth=.5)
ax.set_xlim(xmin, xmax)
ax.set_ylim(ymin, ymax)
ax.add_artist(ScaleBar(1, location="lower left"))
x, y, arrow_length = 0.05, .98, 0.05
ax.annotate('N', xy=(x, y), xytext=(x, y-arrow_length),
arrowprops=dict(facecolor='black', width=5, headwidth=15),
ha='center', va='center', fontsize=20,
xycoords=ax.transAxes)
ax.set_xticks([])
ax.set_yticks([])
1699039.1745048524 1809972.3879951476 5864874.988694763 6008322.51060524
[]